wip

fluid/fluid_problem.c

@@ -122,7 +122,6 @@ void fluid_problem_node_force_volume(FluidProblem *problem, const double *soluti
     }
     double f[D],dfdu,dfddp;
     particle_force_f(problem,f,&dfdu,&dfddp,s,ds,sold,c,dc,a,dt,iel,ip);
-    if(ip == 0) printf("Force : %14.7f\n", f[0]);
 
     if (!local_vector)
       continue;
@@ -372,6 +371,15 @@ static void particle_force_f(FluidProblem *problem, double *f, double *dfdu, dou
   }
   *dfdu = gammastar/c;
   *dfddp = gammastar*dt/mass*vol;//-vol;
+  if(ip == 2){
+    for(int i = 0; i < D; ++i){
+      printf("u : %14.7f\n", u[i]);
+    }
+    printf("c : %14.7f\n", c);
+    printf("f : %14.7f\n", f[0]);
+    printf("dfdu : %14.7f\n", *dfdu);
+    printf("dfddp : %14.7f\n", *dfddp);
+  }
 }
 
 static void fluid_problem_f(const FluidProblem *problem, const double *sol, double dsol[][D], const double *solold, const double dsolold[][D], const double* mesh_velocity, const double c, const double *dc, const double *bf, const double cold, const double rho, const double mu, double dt, int iel, double *f0, double f1[][D], double f00[], double f10[][D], double f01[][D], double f11[][D][D]) {