Merge Gizmo-SPH into the master branch
Merge request reports
Activity
@bvandenbroucke, can I remove
rhofrom the particle definition ? How aboutdu_dt? Do you useomegaandv_fullin thexpart?
omegais not used, so that can go.v_fullis used indirectly inhydro_iact(meaning it is copied into a similar field inpart, sincexpartis not available insidehydro_iact).du_dtis used for the drift. You are welcome to removerho, if the code still compiles if you do that (which it probably won't, sincerhois used in the common part of the drift).Edited by Bert VandenbrouckeAdded 1 commit:
- 8db08c24 - Removed unnecessary variables in the GIZMO particle definition.
Added 1 commit:
- a0cc3bf9 - Common interface to access the mass and the density of particles.
Added 1 commit:
- a3d0c6bc - Also update the prediction set interface fot Gizmo.
Added 1 commit:
- 85612470 - Also update the prediction set interface fot Gizmo.
Added 1 commit:
- 72945cae - Make sure the fluid quantities are converted correctly from ICs to internal stuff in GIZMO.
-
@bvandenbroucke, one last thing, could you check that you are happy with the changes in the last commit 72945cae ? I had to change the energy variable read in as I had removed the mass field. It all runs but just want to make sure you are fine with it.
254 253 __attribute__((always_inline)) INLINE static void hydro_convert_quantities( 255 254 struct part* p) { 256 255 257 float momentum[3]; 258 256 const float volume = p->geometry.volume; 259 257 const float m = p->conserved.mass; 260 258 p->primitives.rho = m / volume; 261 259 262 /* P actually contains internal energy at this point */ 263 p->primitives.P *= hydro_gamma_minus_one * p->primitives.rho; 260 p->conserved.momentum[0] = m * p->primitives.v[0]; 261 p->conserved.momentum[1] = m * p->primitives.v[1]; 262 p->conserved.momentum[2] = m * p->primitives.v[2]; 263 p->conserved.energy *= m; 254 253 __attribute__((always_inline)) INLINE static void hydro_convert_quantities( 255 254 struct part* p) { 256 255 257 float momentum[3]; 258 256 const float volume = p->geometry.volume; 259 257 const float m = p->conserved.mass; 260 258 p->primitives.rho = m / volume; 261 259 262 /* P actually contains internal energy at this point */ 263 p->primitives.P *= hydro_gamma_minus_one * p->primitives.rho; 260 p->conserved.momentum[0] = m * p->primitives.v[0]; 261 p->conserved.momentum[1] = m * p->primitives.v[1]; 262 p->conserved.momentum[2] = m * p->primitives.v[2]; 263 p->conserved.energy *= m; 254 253 __attribute__((always_inline)) INLINE static void hydro_convert_quantities( 255 254 struct part* p) { 256 255 257 float momentum[3]; 258 256 const float volume = p->geometry.volume; 259 257 const float m = p->conserved.mass; 260 258 p->primitives.rho = m / volume; 261 259 262 /* P actually contains internal energy at this point */ 263 p->primitives.P *= hydro_gamma_minus_one * p->primitives.rho; 260 p->conserved.momentum[0] = m * p->primitives.v[0]; 261 p->conserved.momentum[1] = m * p->primitives.v[1]; 262 p->conserved.momentum[2] = m * p->primitives.v[2]; 263 p->conserved.energy *= m; -
No, in the last version that I have the internal energy (not per unit mass) is calculated from the pressure. The pressure initially contains the internal energy per unit mass and is then converted into the pressure.
If u is the internal energy per unit mass, then P = (gamma-1) * rhou; conserved.energy = P / (gamma-1) * volume = um.
What you have now is conserved.energy = um (which is correct), and P = (gamma-1)rhomu (which is not correct). If you put the p->conserved.energy *= m below the computation of P, it will be correct.
Edited by Bert Vandenbroucke
254 253 __attribute__((always_inline)) INLINE static void hydro_convert_quantities( 255 254 struct part* p) { 256 255 257 float momentum[3]; 258 256 const float volume = p->geometry.volume; 259 257 const float m = p->conserved.mass; 260 258 p->primitives.rho = m / volume; 261 259 262 /* P actually contains internal energy at this point */ 263 p->primitives.P *= hydro_gamma_minus_one * p->primitives.rho; 260 p->conserved.momentum[0] = m * p->primitives.v[0]; 261 p->conserved.momentum[1] = m * p->primitives.v[1]; 262 p->conserved.momentum[2] = m * p->primitives.v[2]; 263 p->conserved.energy *= m;