From: David Knezevic <dknezevic@se...>  20130401 13:06:09

On 04/01/2013 07:50 AM, Kirk, Benjamin (JSCEG311) wrote: > On Mar 31, 2013, at 11:35 PM, "Manav Bhatia" <bhatiamanav@...> wrote: > >> Thanks, Derek. >> >> A related question about libMesh: >> >> If I need to output the solution values on the boundary nodes, can I use BoundaryMesh along with the IO classes for that? >> >> Meaning, if I do the following >> >> BoundaryMesh b_mesh(mesh.mesh_dimension()1); >> mesh.boundary_info>sync(b_mesh); >> VTKIO(b_mesh).write_equation_systems("boundary_output.pvtu", equation_systems); >> >> >> Would this write the solution values on the boundary nodes? > Almost. Unfortunately it is a little more manual than this, but that's the right idea. In your example the equation systems object can only be associated with one mesh, so it won't work. The fix is to create a second equation systems object that lives on the boundary and then extract the trace values. > > This is not automated because in general the user will want physicsdependent things like shear stress, heat flux, whatever... > > I suppose we could automate the valuesonly case though. > > Ben To get the "boundary only values" I first make a "boundary_dofmap", as in the code below. Then copying the solution values over from system to boundary_system is easy. void create_disp_boundary_dofmap(System& system, System& boundary_system, std::vector<unsigned int>& boundary_dofmap) { boundary_dofmap.resize(boundary_system.n_dofs()); // make a point locator for system AutoPtr<PointLocatorBase> point_locator = system.get_mesh().sub_point_locator(); // loop over the boundary nodes and locate them // then fill in boundary_dofmap MeshBase::node_iterator node_it = boundary_system.get_mesh().nodes_begin(); const MeshBase::node_iterator node_end = boundary_system.get_mesh().nodes_end(); for( ; node_it != node_end; node_it++) { Node* node = *node_it; // get an element in the full mesh that contains node const Elem* element = point_locator>operator()(*node); // loop over the nodes of element until we find the one that matches node for(unsigned int node_id=0; node_id<element>n_nodes(); node_id++) { Node* new_node = element>get_node(node_id); Real dist = std::sqrt( pow(node>operator()(0)new_node>operator()(0),2.) + pow(node>operator()(1)new_node>operator()(1),2.) + pow(node>operator()(2)new_node>operator()(2),2.) ); if(dist < TOLERANCE) { for(unsigned int var=0; var<system.n_vars(); var++) { unsigned int index1 = node>dof_number(boundary_system.number(), var, 0); unsigned int index2 = new_node>dof_number(system.number(), var, 0); boundary_dofmap[index1] = index2; } continue; } } } } 