init porosity

src/fluid_problem.c

@@ -480,6 +480,11 @@ FluidProblem *fluid_problem_new(const char *mesh_file_name, double g, double mu,
   problem->epsilon = malloc(mesh->n_nodes*sizeof(double));
   problem->old_porosity = malloc(mesh->n_nodes*sizeof(double));
   problem->solution = malloc(mesh->n_nodes*n_fields*sizeof(double));
+  for (int i = 0; i < problem->mesh->n_nodes; ++i)
+  {
+    problem->porosity[i] = 1.;
+    problem->old_porosity[i] = 1.;
+  }
   // begin to remove
   for (int i = 0; i < problem->mesh->n_nodes*N_N; ++i){
     problem->solution[i] = 0.;
@@ -763,6 +768,11 @@ void fluid_problem_adapt_mesh(FluidProblem *problem, double lcmax, double lcmin,
   problem->porosity = malloc(sizeof(double)*new_mesh->n_nodes);
   free(problem->old_porosity);
   problem->old_porosity = malloc(sizeof(double)*new_mesh->n_nodes);
+  for (int i = 0; i < new_mesh->n_nodes; ++i)
+  {
+    problem->porosity[i] = 1.;
+    problem->old_porosity[i] = 1.;
+  }
   mesh_free(problem->mesh);
   problem->mesh = new_mesh;
   for (int i = 0; i < problem->n_particles; ++i)