#6 Bug in initializing velocity field
I find the bug to explain why we have extremely large velocity at the bottom cell.
In the example you sent to me, we initialize the velocity field with prop->n1. However, if you check the InitializePhysicalVariables function in phys.c, you will find that we use grid->Nkc not grid->Nke to set initial condition of velocity field. Regularly, it's not a problem, however when Nke<Nkc, we will have a big problem.
Here is the example. If Nke=8 and Nkc=9, because we initialize our velocity field with Nkc, so we will have a fake Initial velocity at u[ne][8]. And because Nke=8, so the velocity of u[ne][8] never changes and keeps as initial velocity. This fake velocity has nothing to do with mass conservation due to the zero flux height. However, when we calculate uc[i][8] according to Nkc (see function ComputeUCperot), we will use the velocity u[ne][8] which should be zero. The non-zero constant u[ne][8] will make a fake uc results based on the initial condition.
In order to solve this bug, you can simply change the way to set IC for u from 0 to Nke[ne] in InitializePhysicalVariables (phys.c line 624-633).
// Initialize the velocity field
for(j=0;j<grid->Ne;j++) {
z = 0;
for(k=0;k<grid->Nke[j] <change from="" nkc<span="">[j]> ;k++) {
z-=grid->dz[k]/2;
phys->u[j][k]=ReturnHorizontalVelocity(
grid->xe[j],grid->ye[j],grid->n1[j],grid->n2[j],z);
z-=grid->dz[k]/2;
}
}</change></grid-></grid->
