[mpas-developers] Vector Reconstruction

Todd Ringler ringler at lanl.gov
Wed Mar 24 23:20:09 MDT 2010


Hi Michael,

I would be in favor or rewriting the RBF bit of code and addressing  
assumptions in the grid structures at a later point.

Cheers,
Todd


On Mar 24, 2010, at 4:30 PM, Michael Duda wrote:

> Hi, Todd.
>
> You're right: removing -D_MPI should work to run the code in
> serial. I'd agree that moving this to a definition at the top of
> the Makefile would be a good way to simplify the work of compiling
> the code for serial execution.
>
> Regarding the differences in behavior of the vector reconstruction
> when changing the grid.nc file, I think I may have found the
> issue. The vector reconstruction code makes the assumption that
> edgesOnCell(i,iCell) is the edge separating iCell from
> cellsOnCell(i,iCell), i=1,nEdgesOnCell(iCell), and this assumption
> doesn't appear to hold for the grid.nc file that I've been testing
> with (based on some debugging code I added to the vector
> reconstruction module). One solution would be to rewrite the code
>
>         do j=1,nEdgesOnCell(iCell)
>           cell2 = cellsOnCell(j,iCell)
>           iEdge = edgesOnCell(j,iCell)
>           X2(1) = xCell(cell2)
>           X2(2) = yCell(cell2)
>           X2(3) = zCell(cell2)
>           if (iCell == cellsOnEdge(1,iEdge)) then
>               normals(:,j) = X2(:) - X1(:)
>           else
>               normals(:,j) = X1(:) - X2(:)
>           endif
>
>           call unit_vector_in_R3(normals(:,j))
>
>           xLoc(:,j,iCell) = 0.5*(X2(:) + X1(:))
>         enddo
>
> to
>
>         do j=1,nEdgesOnCell(iCell)
>           iEdge = edgesOnCell(j,iCell)
>           if (iCell == cellsOnEdge(1,iEdge)) then
>               cell2 = cellsOnEdge(2,iEdge)
>               X2(1) = xCell(cell2)
>               X2(2) = yCell(cell2)
>               X2(3) = zCell(cell2)
>               normals(:,j) = X2(:) - X1(:)
>           else
>               cell2 = cellsOnEdge(1,iEdge)
>               X2(1) = xCell(cell2)
>               X2(2) = yCell(cell2)
>               X2(3) = zCell(cell2)
>               normals(:,j) = X1(:) - X2(:)
>           endif
>
>           call unit_vector_in_R3(normals(:,j))
>
>           xLoc(:,j,iCell) = 0.5*(X2(:) + X1(:))
>         enddo
>
>
> to remove the assumption; doing this, I get reasonable values for
> uReconstructZ for SWTC2. Another solution might be to just require
> that grid.nc files provide edgesOnCell and cellsOnCell arrays that
> begin at the same angle and list edges/cells in counter-clockwise
> order, satisfying the assumption currently made in the vector
> reconstruction code; it sounds like our current repository
> versions of the grid generators might already satisfy this, but it
> would be something that we'd need to specify in case grids were to
> be generated through alternative codes. Given that the fix in the
> vector reconstruction code is relatively small, I'd lean in favor
> of imposing fewer requirements on grid.nc files at this point, but
> we might find that, down the road a bit, having the luxury of
> making this assumption would make our lives much easier. What do
> you all think?
>
> Cheers!
> Michael
>
>
> On Wed, Mar 24, 2010 at 04:00:43PM -0600, Todd Ringler wrote:
>>
>> Michael,
>>
>> When I want to run in serial mode I just remove the "-D_MPI" compile
>> flag under the trunk code. I am not sure if this is suppose to work,
>> but it does (at least for the times I have tried it).
>>
>> I recompiled the trunk sw for serial model and I get errors in z
>> velocity of about 8e-3.
>>
>> (If we think that the trunk code should work in serial mode when  
>> D_MPI
>> is omitted, maybe at some point we can move that option to the top of
>> the Makefile to allow folks to commit it out.)
>>
>> Cheers,
>> Todd
>>
>>
>> On Mar 24, 2010, at 3:33 PM, Michael Duda wrote:
>>
>>> Hi, Todd.
>>>
>>> This is definitely very odd. The executable is named swmodel
>>> rather than sw_model.exe since the code is from the swmodel branch
>>> (formerly in the trunk), which doesn't support multiple cores
>>> under the same directory structure. In our current trunk, the
>>> executable is always named $(CORE)_model.exe, whereas we had
>>> simply called it swmodel in the previous trunk/swmodel code.
>>>
>>> The grid is the quasi-uniform 10242 mesh that I've had posted at
>>> http://www.mmm.ucar.edu/people/duda/files/mpas, and was created
>>> based on the generating points you had sent sometime late last
>>> year. I don't think the problem is connected with the mesh, since
>>> the 'u' and 'v' fields (i.e., the normal and tangential velocity
>>> fields) in the output.nc file look correct when I plot them with
>>> ncl/cells.ncl. The height field from the simulation remains
>>> stable, so I don't think there's a problem with the underlying
>>> simulation.
>>>
>>> If you or Xylar could place a tar file of the code, either
>>> trunk/mpas or branches/swmodel, that you've used in your tests
>>> with the new module_vector_reconstruction.F, I could peruse that
>>> and see whether I can spot any key differences in how the module
>>> is used. Just let me know if there's any testing that I can do to
>>> help out here.
>>>
>>> Cheers,
>>> Michael
>>>
>>>
>>> On Wed, Mar 24, 2010 at 03:24:58PM -0600, Todd Ringler wrote:
>>>>
>>>> Hi Michael,
>>>>
>>>> I also reproduced your results.
>>>>
>>>> How can this be? We must be missing something obvious here.
>>>>
>>>> I notice that the executable is named "swmodel" and sw_model.exe --
>>>> was the reconstruction tested on top of the trunk? Also, is this  
>>>> the
>>>> quasi uniform 10242 mesh? built with the current version of the  
>>>> grid
>>>> generator?
>>>>
>>>> Cheers,
>>>> Todd
>>>>
>>>> On Mar 24, 2010, at 3:00 PM, Xylar S. Asay-Davis wrote:
>>>>
>>>>> Hi Michael,
>>>>>
>>>>> I'm able to reproduce your results.  It doesn't look good!  I'll  
>>>>> see
>>>>> what
>>>>> I can figure out.
>>>>>
>>>>> -Xylar
>>>>>
>>>>>
>>>>>> Hi, Todd.
>>>>>>
>>>>>> I wonder whether there's something odd with how I've added the  
>>>>>> new
>>>>>> module_vector_reconstruction.F code and Registry entry. I had  
>>>>>> been
>>>>>> testing with branches/swmodel for simplicity, and the changes are
>>>>>> really trivial. If you wouldn't mind, could you give a try with
>>>>>> the tar file I've been testing with? I've placed this at
>>>>>> http://www.mmm.ucar.edu/people/duda/files/swmodel_vectors.tar.gz.
>>>>>> The namelist is set up to run just 5 time steps, and the input
>>>>>> files are all included, so it should just be a matter of  
>>>>>> compiling
>>>>>> and running. Also included in the tar file is the output file  
>>>>>> that
>>>>>> I generated for SWTC2; doing an ncdump or using ncview, you can
>>>>>> see the large values in uReconstructZ.
>>>>>>
>>>>>> If you're able to get correct results with my tar file, I'd not
>>>>>> object to committing the new module_vector_reconstruction.F; I  
>>>>>> can
>>>>>> propose to add this in a new shared directory 'operators', to be
>>>>>> created in trunk/mpas/src with appropriate changes to the
>>>>>> Makefiles. However, if you can reproduce the large values of
>>>>>> uReconstructZ with my tar file, I can look for any differences if
>>>>>> you could supply the code that you were able to successfully run
>>>>>> with. Does that sound like a reasonable way forward?
>>>>>>
>>>>>> Cheers,
>>>>>> Michael
>>>>>>
>>>>>>
>>>>>> On Wed, Mar 24, 2010 at 12:06:19PM -0600, Todd Ringler wrote:
>>>>>>>
>>>>>>> Hi Michael,
>>>>>>>
>>>>>>> I can not reproduce your finding. Below is the rbf  
>>>>>>> reconstruction
>>>>>>> (vectors) and z-component of velocity (colors with colorbar,  
>>>>>>> units
>>>>>>> m/
>>>>>>> s). This is SWTC#2 on the quasi-uniform 2562 grid. I am happy to
>>>>>>> followup over the phone. Have to run, back at 2 pm.
>>>>>>>
>>>>>>> Cheers,
>>>>>>> Todd
>>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>>>
>>>>>>> On Mar 23, 2010, at 11:56 AM, Michael Duda wrote:
>>>>>>>
>>>>>>>> Hi, Xylar.
>>>>>>>>
>>>>>>>> Thanks for the quick reply. In SWTC2, the maximum zonal  
>>>>>>>> velocity
>>>>>>>> is around 38 m/s at the equator, so 10 m/s for uReconstructZ  
>>>>>>>> does
>>>>>>>> seem out of place. Just to make sure I understand how to
>>>>>>>> interpret
>>>>>>>> uReconstruct{X,Y,Z}, these are the components of the velocity
>>>>>>>> vector in Cartesian space, with the z-axis positive through the
>>>>>>>> north pole of the sphere, right? I had been testing on a fairly
>>>>>>>> coarse grid (10242 cells), so that could be a factor; I'll do a
>>>>>>>> little more testing at higher resolution.
>>>>>>>>
>>>>>>>> Cheers,
>>>>>>>> Michael
>>>>>>>>
>>>>>>>>
>>>>>>>> On Tue, Mar 23, 2010 at 11:41:18AM -0600, Xylar S. Asay-Davis
>>>>>>>> wrote:
>>>>>>>>> Hi Michael,
>>>>>>>>>
>>>>>>>>> I don't think your attachment made it through.
>>>>>>>>>
>>>>>>>>> I'll take a look at shallow water test case 2 to see what
>>>>>>>>> might be
>>>>>>>>> the
>>>>>>>>> problem there. The value of uReconstructZ won't machine zero
>>>>>>>>> for a
>>>>>>>>> purely
>>>>>>>>> zonal flow because the reconstruction in its current form will
>>>>>>>>> only
>>>>>>>>> guarantee exact reconstruct of a constant flow in the tangent
>>>>>>>>> plane
>>>>>>>>> of the
>>>>>>>>> cell, where as a zonal flow has components slightly out of the
>>>>>>>>> tangent
>>>>>>>>> plane.  In my tests, I was seeing uReconstructZ smaller than
>>>>>>>>> uReconstructX
>>>>>>>>> and uReconstructY by a factor of around 1e-6.
>>>>>>>>>
>>>>>>>>> -Xylar
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>> Hi, Xylar.
>>>>>>>>>>
>>>>>>>>>> The module looks good to me -- very clean code. The only
>>>>>>>>>> changes
>>>>>>>>>> I'd suggest are to the comments at the beginning of the
>>>>>>>>>> init_reconstruct() routine, and a few changes to indentation,
>>>>>>>>>> as
>>>>>>>>>> in the attached module_vector_reconstruction.F.
>>>>>>>>>>
>>>>>>>>>> This is most likely a problem in how I'm performing the
>>>>>>>>>> computation given that both you and Todd have checked the
>>>>>>>>>> reconstructed vectors, but I'm having trouble getting a zonal
>>>>>>>>>> and
>>>>>>>>>> meridional wind from the reconstructed vector in (x,y,z).  
>>>>>>>>>> Even
>>>>>>>>>> plotting the magnitude of the reconstructed vector,  
>>>>>>>>>> sqrt(x^2 +
>>>>>>>>>> y^2
>>>>>>>>>> + z^2), I'm getting some unexpected values in the  
>>>>>>>>>> reconstructed
>>>>>>>>>> wind speed field for shallow water test case 2 (with  
>>>>>>>>>> alpha=0).
>>>>>>>>>> Given that we should have purely zonal flow for this test  
>>>>>>>>>> case,
>>>>>>>>>> I'm suprised to see values in uReconstructZ of order 10. Do
>>>>>>>>>> either
>>>>>>>>>> you or Todd have any thoughts or suggestions?
>>>>>>>>>>
>>>>>>>>>> Cheers,
>>>>>>>>>> Michael
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> On Mon, Mar 22, 2010 at 03:30:39PM -0600, Xylar S. Asay-Davis
>>>>>>>>>> wrote:
>>>>>>>>>>> Attached is my proposed vector reconstruction code.  I have
>>>>>>>>>>> run
>>>>>>>>>>> several
>>>>>>>>>>> sanity checks and it seems to work fine.  However, I haven't
>>>>>>>>>>> been
>>>>>>>>>>> able
>>>>>>>>>>> to
>>>>>>>>>>> get a full ocean core run with this code.  It is not clear
>>>>>>>>>>> to me
>>>>>>>>>>> if this
>>>>>>>>>>> has to do with this code or with recent changes that Mark  
>>>>>>>>>>> has
>>>>>>>>>>> made to
>>>>>>>>>>> the
>>>>>>>>>>> ocean core relating to 3D density.  I'll keep working with
>>>>>>>>>>> Mark
>>>>>>>>>>> to try
>>>>>>>>>>> to
>>>>>>>>>>> figure that out, but in the mean time it would be great to
>>>>>>>>>>> have
>>>>>>>>>>> Todd
>>>>>>>>>>> and/or Michael test out this code to see if it will work for
>>>>>>>>>>> you.
>>>>>>>>>>>
>>>>>>>>>>> A registry entry for the coeffs_reconst is required as  
>>>>>>>>>>> follows
>>>>>>>>>>> (already
>>>>>>>>>>> in
>>>>>>>>>>> the ocean core):
>>>>>>>>>>> var real    coeffs_reconstruct ( R3 maxEdges nCells ) -
>>>>>>>>>>> coeffs_reconstruct
>>>>>>>>>>> - -
>>>>>>>>>>>
>>>>>>>>>>> Registry entries for  matrix_reconstruct, normal and
>>>>>>>>>>> rbf_value
>>>>>>>>>>> are no
>>>>>>>>>>> longer required by this code.
>>>>>>>>>>>
>>>>>>>>>>> -Xylar
>>>>>>>>>>
>>>>>>>>>>> _______________________________________________
>>>>>>>>>>> mpas-developers mailing list
>>>>>>>>>>> mpas-developers at mailman.ucar.edu
>>>>>>>>>>> http://mailman.ucar.edu/mailman/listinfo/mpas-developers
>>>>>>>>>>
>>>>>>>>>> _______________________________________________
>>>>>>>>>> mpas-developers mailing list
>>>>>>>>>> mpas-developers at mailman.ucar.edu
>>>>>>>>>> http://mailman.ucar.edu/mailman/listinfo/mpas-developers
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> ***********************
>>>>>>>>> Xylar S. Asay-Davis
>>>>>>>>> E-mail: xylar at lanl.gov
>>>>>>>>> Phone: (505) 606-0025
>>>>>>>>> Fax: (505) 665-2659
>>>>>>>>> CNLS, MS B258
>>>>>>>>> Los Alamos National Laboratory
>>>>>>>>> Los Alamos, NM 87545
>>>>>>>>> ***********************
>>>>>>>>>
>>>>>>>> _______________________________________________
>>>>>>>> mpas-developers mailing list
>>>>>>>> mpas-developers at mailman.ucar.edu
>>>>>>>> http://mailman.ucar.edu/mailman/listinfo/mpas-developers
>>>>>>>
>>>>>>
>>>>>> _______________________________________________
>>>>>> mpas-developers mailing list
>>>>>> mpas-developers at mailman.ucar.edu
>>>>>> http://mailman.ucar.edu/mailman/listinfo/mpas-developers
>>>>>>
>>>>>
>>>>>
>>>>> ***********************
>>>>> Xylar S. Asay-Davis
>>>>> E-mail: xylar at lanl.gov
>>>>> Phone: (505) 606-0025
>>>>> Fax: (505) 665-2659
>>>>> CNLS, MS B258
>>>>> Los Alamos National Laboratory
>>>>> Los Alamos, NM 87545
>>>>> ***********************
>>>>>
>>>>
>>
> _______________________________________________
> mpas-developers mailing list
> mpas-developers at mailman.ucar.edu
> http://mailman.ucar.edu/mailman/listinfo/mpas-developers



More information about the mpas-developers mailing list