update_tracer.cpp

Go to the documentation of this file.

#include "iMesh.h"
#include "MBiMesh.hpp"
#include "moab/Core.hpp"
#include "moab/Range.hpp"
#include "moab/IntxMesh/Intx2MeshOnSphere.hpp"
#include "moab/IntxMesh/IntxUtils.hpp"
#include "IntxUtilsCSLAM.hpp"
 
extern "C" void update_tracer( iMesh_Instance instance, iBase_EntitySetHandle imesh_euler_set, int* ierr )
{
 using namespace moab;
 const double radius = 1.;
 const double gtol = 1.e-9;
 const bool debug = false;
 
 Range ents;
 moab::Interface* mb = MOABI;
 *ierr = 1;
 
 EntityHandle euler_set = (EntityHandle)imesh_euler_set;
 
 Intx2MeshOnSphere worker( mb );
 worker.set_radius_source_mesh( radius );
 worker.set_radius_destination_mesh( radius );
 worker.set_error_tolerance( gtol );
 
 EntityHandle covering_lagr_set;
 
 ErrorCode rval = mb->create_meshset( MESHSET_SET, covering_lagr_set );MB_CHK_SET_ERR_RET( rval, "can't create covering set " );
 
 // we need to update the correlation tag and remote tuples
 rval = worker.create_departure_mesh_2nd_alg( euler_set, covering_lagr_set );MB_CHK_SET_ERR_RET( rval, "can't populate covering set " );
 
 if( debug )
 {
 rval = mb->write_file( "lagr.h5m", 0, 0, &covering_lagr_set, 1 );MB_CHK_SET_ERR_RET( rval, "can't write covering set " );
 }
 
 // enforce convexity for the set
 rval = moab::IntxUtils::enforce_convexity( mb, covering_lagr_set );MB_CHK_SET_ERR_RET( rval, "can't write covering set " );
 
 EntityHandle outputSet;
 rval = mb->create_meshset( MESHSET_SET, outputSet );MB_CHK_SET_ERR_RET( rval, "can't create output set " );
 
 rval = worker.intersect_meshes( covering_lagr_set, euler_set, outputSet );MB_CHK_SET_ERR_RET( rval, "can't intersect " );
 
 if( debug )
 {
 rval = mb->write_file( "output.vtk", 0, 0, &outputSet, 1 );MB_CHK_SET_ERR_RET( rval, "can't write covering set " );
 }
 
 // tagElem is the average computed at each element, from nodal values
 Tag tagElem = 0;
 std::string tag_name2( "TracerAverage" );
 rval = mb->tag_get_handle( tag_name2.c_str(), 1, MB_TYPE_DOUBLE, tagElem, MB_TAG_DENSE | MB_TAG_CREAT );MB_CHK_SET_ERR_RET( rval, "can't get tracer tag " );
 
 // area of the euler element is fixed, store it; it is used to recompute the averages at each
 // time step
 Tag tagArea = 0;
 std::string tag_name4( "Area" );
 rval = mb->tag_get_handle( tag_name4.c_str(), 1, MB_TYPE_DOUBLE, tagArea, MB_TAG_DENSE | MB_TAG_CREAT );MB_CHK_SET_ERR_RET( rval, "can't get area tag" );
 
 rval = worker.update_tracer_data( outputSet, tagElem, tagArea );MB_CHK_SET_ERR_RET( rval, "can't update tracer " );
 
 // everything can be deleted now from intx data; polygons, etc.
 
 *ierr = 0;
 return;
}