old stab

fluid/fluid_problem.c

@@ -565,7 +565,7 @@ static void fluid_problem_f(FluidProblem *problem, double *f0, double *f1, doubl
       f1[(U+i)*D+j] = -rho*u[i]*u[j]/c + mu*(tau[i][j]+tau[j][i]) + (i==j ? -p : 0) + R*u[j]*rho*taup;
     }
     f1[(U+i)*D+i] += divu*tauc*rho;
-    f1[P*D+i] = -u[i] + 5*taup*R;
+    f1[P*D+i] = -u[i] + 0.025*dp[i];
     if(problem->n_fluids == 2){
       f1[P*D+i] +=  taua*(divu+(c-cold)/dt)*u[i];
       f1[A*D+i] = -u[i]*a + taua*Ra*u[i] + taup*R*a + da[i]*(1e-7+extraa);
@@ -883,7 +883,7 @@ static void fluid_problem_assemble_system(FluidProblem *problem, double *rhs, co
     }
   }
   node_force_volume(problem, solution_old, dt, all_local_vector, all_local_matrix,reduced_gravity);
-  //compute_weak_boundary_conditions(problem, dt, all_local_vector, all_local_matrix);
+  compute_weak_boundary_conditions(problem, dt, all_local_vector, all_local_matrix);
   fluid_problem_volume(problem, solution_old, dt, all_local_vector, all_local_matrix,reduced_gravity);
   for (int iel=0; iel < mesh->n_elements; ++iel) {
     double *local_vector = &all_local_vector[local_size*iel];