<html><body style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space;" class="">Hello,<br class=""><br class="">I am learning to use Petsc in order to modify my sequential FEM code to work in parallel. First, I am trying to understand how a FEM assembly loop would work in parallel using Petsc but I run into some problems. <br class=""><br class=""><div class="">I made a simple example to explain the issue : <br class="">- I read a simple mesh (generated by gmsh) with 4 quadrilaterls and generate the corresponding DMPlex<br class="">- I define a scalar field on vertices <br class="">- Vertices on Bottom and Top boundaries are contrained (so I have 3 unknowns)<br class="">- I create a DM Matrix <br class="">- For the purpose of my example during the assembly I define an “elementary matrix” containing only ones <br class=""><br class="">When I run the program with one process I obtain what the matrix is supposed to be : <br class=""><br class="">Mat Object: 1 MPI processes<br class=""> type: seqaij<br class="">row 0: (0, 2.) (2, 2.)<br class="">row 1: (1, 2.) (2, 2.)<br class="">row 2: (0, 2.) (1, 2.) (2, 4.)<br class=""><br class="">But running with two processes gives me this :<br class=""><br class="">Mat Object: 2 MPI processes<br class=""> type: mpiaij<br class="">row 0: (0, 2.) (2, 0.)<br class="">row 1: (1, 2.) (2, 2.)<br class="">row 2: (0, 0.) (1, 2.) (2, 2.)<br class=""><br class="">It seems that in the way I am doing the assembly, Process 0 is not giving the contribution of its local assembly to process 1. <div class=""><br class=""></div><div class=""><br class="">For the assembly I am doing something like this (the complete example is attached) : <br class=""><br class=""><br class=""><i class=""> !*************************************<br class=""> !--- Matrix Assembly<br class=""> !*************************************<br class=""><br class=""> call DMCreateMatrix(dm,A,ierr)<br class=""><span class=""> call DMGetLocalToGlobalMapping(dm,ltog,ierr)<br class=""></span><span class=""> call MatSetLocalToGlobalMapping(A,ltog,ltog,ierr);<br class=""></span><span class=""> call ISLocalToGlobalMappingDestroy(ltog,ierr);<br class=""></span><br class=""><span class=""><br class=""></span><span class="">do icell = </span>!--- Loop over cells<span class=""><br class=""></span><span class=""> </span><span class=""><br class=""> do ic = ! Loop over closure (somehow only vertices)</span><span class=""><br class=""> ipoin = closure((ic-1)*2+1)<br class=""></span><span class=""> call DMPlexGetPointLocal(dm,ipoin,pStart,pEnd,ierr)<br class=""> ! <br class=""> ! Something to get location of point data in a local Vec <br class=""> ! <br class=""> enddo<br class=""><br class=""> elstf(:,:) = 1.0<br class=""> </span></i><div class=""><span class=""><i class=""> do i = !Loop over rows corresponding to dofs which are not constrained<br class=""> call MatSetValuesLocal(A,1,local_ldof(i)-1,ndof*nodes,local_ldof(:)-1, elstf(i,:),ADD_VALUES,ierr)<br class=""> enddo</i></span></div><div class=""><span class=""><i class=""> <br class="">enddo<br class=""> <br class=""> call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)<br class=""> call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)</i></span></div><div class=""><span class=""><i class=""><br class=""></i></span></div><div class=""><i class=""> !*************************************<br class=""> !--- End Assembly<br class=""> !*************************************</i></div><div class=""><span class=""><br class=""></span></div><div class=""><span class=""><br class=""></span></div><div class=""><br class=""></div><div class=""><span class="">There is something I am definitely missing because I don’t understand why it does not work. I was thinking maybe I need a function like </span><span class="">DMLocallToGlobal…() but for matrices. </span></div><span class="">I hope someone can help me with this and I am sorry if there is an obvious answer</span>…</div><div class=""><span class="">Thanks in advance. </span></div><div class=""><span class=""><br class=""></span></div><div class=""><span class="">Nestor<br class=""></span><div class=""><span class=""><br class=""></span></div><div class=""><span class=""><br class=""><br class=""></span></div></div></div></body></html>