Thanks. It’s really helpful for me. I have another question about local nodes 
and local elements.

A process can only access one part of solution vectors. This part of vectors 
contains the solution of all local nodes, is it right?
If it’s right, when there is a local element contains a node that is not owned 
by this process, how to get the solution vector of this node without requesting 
other process? 
And can point_value() be executed to get the value of a point located in this 
element since some nodes are not owned by the process? 


> On Jun 27, 2016, at 1:33 PM, John Peterson <jwpeter...@gmail.com> wrote:
> 
> 
> 
> On Mon, Jun 27, 2016 at 10:24 AM, 张江 <zhangjiang.d...@gmail.com 
> <mailto:zhangjiang.d...@gmail.com>> wrote:
> Hi,
> 
> I have a few questions about reading data in parallel using libmesh. The 
> program is running on n processes, and I wish to have n partitions on this 
> mesh which means each data partition is assigned to one process.
> 
> 1. How to make sure that each process just contains one part of the entire 
> mesh and the corresponding nodal values?
> 2. How can each process maintain the data values of ghost nodes since its 
> local elements may contains the nodes that are not in its local nodes?
> 
> Here is the initialization code. It seems that each process can access the 
> entire mesh and the entire nodal solution. Does it mean that each process 
> hold the entire data?
> 
> I'm not sure what reason you have for thinking this... if you look at the 
> code for ExodusII_IO::copy_nodal_solution() you can see that it reads the 
> nodal values on every processor, but it only sets the values for which it is 
> responsible.
> 
> Generally speaking, you should only access dofs which fall in the half-open 
> range
> [system.solution->first_local_index(), system.solution->last_local_index())
> 
>  
> 
> libMesh::UnstructuredMesh *_mesh;
> libMesh::ExodusII_IO *_exio;
> libMesh::EquationSystems *_eqsys;
> libMesh::System *_asys;
> libMesh::PointLocatorTree *_locator;
> 
> _mesh = new Mesh(comm());
> _exio = new ExodusII_IO(*_mesh);
> _exio->read(filename);
> _mesh->allow_renumbering(false);
> _mesh->prepare_for_use();
> 
> General comment: don't use "new" and "delete" unless you have an explicit 
> reason for doing so.
> 
>  
> 
> unsigned int _Vx_var, _Vy_var, _Vz_var;
> _eqsys = new EquationSystems(*_mesh);
> _asys = &(_eqsys->add_system<System>("Velocity"));
> _Vx_var = _asys->add_variable("Vx", FIRST, LAGRANGE);
> _Vy_var = _asys->add_variable("Vy", FIRST, LAGRANGE);
> _Vz_var = _asys->add_variable("Vz", FIRST, LAGRANGE);
> 
> _eqsys->init();
> 
> _locator = new PointLocatorTree(*_mesh);
> 
> _exio->copy_nodal_solution(*_asys, "Vx", varnames[0]); // velocityX
> _exio->copy_nodal_solution(*_asys, "Vy", varnames[1]); // velocityY
> _exio->copy_nodal_solution(*_asys, "Vz", varnames[2]); // velocityZ
> 
> libMesh::UniquePtr<libMesh::NumericVector<libMesh::Number> > _as = 
> _asys->current_local_solution->clone();
> 
> You are making a copy of the entire current_local_solution vector here, which 
> is almost certainly not what you want.
> 
> -- 
> John

------------------------------------------------------------------------------
Attend Shape: An AT&T Tech Expo July 15-16. Meet us at AT&T Park in San
Francisco, CA to explore cutting-edge tech and listen to tech luminaries
present their vision of the future. This family event has something for
everyone, including kids. Get more information and register today.
http://sdm.link/attshape
_______________________________________________
Libmesh-users mailing list
Libmesh-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/libmesh-users

Reply via email to