- Thank you received: 0
extracting solution on a given plane
5 years 7 months ago #2319
by gcdiwan
extracting solution on a given plane was created by gcdiwan
Dear NGSolve Developers,
I am trying to extract the solution in a 3D problem on a given plane. So, once the solution is available from the solver, I set up a cloud of points defined over a 2D region [-R,R]^2:
While this works without mpirun, when i launch the script with more than 1 cores, i end up in SIGSEGV fault. Most likely because inside the for loop the partitioning information is lost. Is there a better way of extracting the solution on a given plane which works with mpirun? Other possibility which i could think of is to write the vtk file and then use vtkcutter to slice through the domain and save the solution.
Thanks in advance!
I am trying to extract the solution in a 3D problem on a given plane. So, once the solution is available from the solver, I set up a cloud of points defined over a 2D region [-R,R]^2:
Code:
comm.Barrier()
if (rank==0):
Nx=100;
Ny=100;
xmin,xmax,ymin,ymax = [-R,R,-R,R]; # R is given
plot_grid = np.mgrid[xmin:xmax:Nx*1j,ymin:ymax:Ny*1j];
points = np.vstack((plot_grid[0].ravel(),plot_grid[1].ravel(),
np.array([0]*plot_grid[0].size)));
fem_xx = my_points[0,:]
fem_xy = my_points[1,:]
for ipt in range(0,len(fem_xx)):
val = u(mesh(fem_xx[ipt],fem_xy[ipt],0.0))
While this works without mpirun, when i launch the script with more than 1 cores, i end up in SIGSEGV fault. Most likely because inside the for loop the partitioning information is lost. Is there a better way of extracting the solution on a given plane which works with mpirun? Other possibility which i could think of is to write the vtk file and then use vtkcutter to slice through the domain and save the solution.
Thanks in advance!
5 years 7 months ago #2320
by lkogler
Replied by lkogler on topic extracting solution on a given plane
With MPI, every proc only has one part of the mesh, with rank 0 not keeping any part at all.
The code
will be called on every rank. Every rank will only have elements for SOME of the points in your plane.
Currently it segfaults for points that are not in the local part of the mesh.
I have attached a small example where I first check every point for whether it is in the mesh or not:
However, in that case you still have your points and values distributed among the MPI procs. If you want them all in one place, mpi4py might be what you are looking for.
Depending on what you want to do with the restriction of your solution, constructing your mesh such that it resolves the plane might be another solution.
Otherwise the VTK approach is also a way.
The code
Code:
for ipt in range(0,len(fem_xx)):
val = u(mesh(fem_xx[ipt],fem_xy[ipt],0.0))
Currently it segfaults for points that are not in the local part of the mesh.
I have attached a small example where I first check every point for whether it is in the mesh or not:
Code:
for ipt in range(0,len(fem_xx)):
if mesh.Contains(fem_xx[ipt],fem_xy[ipt],0.0):
mesh_pnt = mesh(fem_xx[ipt],fem_xy[ipt],0.0)
val = gfu(mesh_pnt)
pvs.append((mesh_pnt, val))
However, in that case you still have your points and values distributed among the MPI procs. If you want them all in one place, mpi4py might be what you are looking for.
Depending on what you want to do with the restriction of your solution, constructing your mesh such that it resolves the plane might be another solution.
Otherwise the VTK approach is also a way.
Attachments:
Time to create page: 0.112 seconds