Preview: Use cell area and volume for finite volume code#3392
Preview: Use cell area and volume for finite volume code#3392dschwoerer wants to merge 35 commits into
Conversation
Information like perpendicular J, Bxyz and g_22 are only known to the grid generator, but are needed for more correct parallel derivatives
_[A-Z].* is preserved for compilers, thus use _[a-z].*
At this time the parallel transform is not yet set, so isFci() return false.
Also specify area for testing
| ASSERT0(mesh->xstart > 0); | ||
| BOUT_FOR(i, area_centre.getRegion("RGN_NOX")) { | ||
| (*_cell_area_xlow)[i] = 0.5 * (area_centre[i] + area_centre[i.xm()]); | ||
| (*_cell_area_xhigh)[i] = 0.5 * (area_centre[i] + area_centre[i.xp()]); |
There was a problem hiding this comment.
The factor 0.5 does imply here that the position of the interface is right at the geometric middle of the cells right? I don't see a reason to not weight this by the actual size of the cell.
There was a problem hiding this comment.
Right now it moves over what the FA code does. I think hypnotoad might have more details about the size and direction of the actual boundary of the cell, but unless hypnotoad adds something like the dagp coefficients, like zoidberg does, it will not be possible to correctly reconstruct it.
So 0.5 is motivated, by several facts
- it is easy
- it is in line with the current FA code
- Anything more complicated immediately needs to go to quite a level of detail, to include extending in the x direction, but also includes shearing, so it needs to depends on g_11, g_13 and g_12
If you have a formulation that is numerically stable, and can show that it gives better results, I am happy to change it. But to me that feels like overkill.
| } | ||
| BOUT_FOR(i, By_c.getRegion("RGN_NOY")) { | ||
| jxz_ylow[i] = By_c[i] / By_l[i] * jxz_centre[i]; | ||
| jxz_yhigh[i] = By_c[i] / By_h[i] * jxz_centre[i]; |
There was a problem hiding this comment.
@dschwoerer We talked about this, again I am not sure that we are taking into account the "proper" conservation of magnetic flux, By here is just the toroidal field and not the total field strength. Maybe @bendudson or @bshanahan have an idea on how to handle this properly. I would argue that when the field line travels into an area of the same toroidal but a different poloidal field, this effect should have an impact on the cell area but here it does not. I know that this is tricky to do because right now we assume that perp = poloidal.
There was a problem hiding this comment.
Lets iterate the same arguments another round :-)
By is the part that is orthogonal to the plane.
If we want to follow the magnetic flux tube, we need to take the area orthogonal to the magnetic field.
That area in the RZ cross section area_{rz} - and the area orthogonal to the magnetic field is smaller by a factor of By / Bmag, a sign of By changes the direction of the area. Thus the total flux would be Bmag * By / Bmag * area_xz = By * area_xz.
I do not see any further assumptions needed for using this expressions.
Certainly, we make further assumptions, but I think if they need to be corrected here in this expressions, we should just remove those assumptions, rather then making it wrong here, to correct for some assumptions some where else.
Note that this is also correct within the assumption that By = Bmag.
| BOUT_FOR(i, area_centre.getRegion("RGN_NOZ")) { | ||
| (*_cell_area_zlow)[i] = 0.5 * (area_centre[i] + area_centre[i.zm()]); | ||
| (*_cell_area_zhigh)[i] = 0.5 * (area_centre[i] + area_centre[i.zp()]); | ||
| } |
0651913 to
3b0117e
Compare
A pr of #3326 into #3335 to show the facility of cell volumes and areas