<div dir="ltr"><div dir="ltr">I think you would prefer to preallocate the matrix<br><br></div><div>uncomment this line <br></div><div dir="ltr"><div>! call MatMPIAIJSetPreallocation(Amat,f9,PETSC_NULL_INTEGER,f6,PETSC_NULL_INTEGER, ierr)<br><br><br></div></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail-m_-5873209276696769746gmail_attr">Le mar. 29 janv. 2019, à 12 h 40, Yaxiong Chen via petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a>> a écrit :<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div dir="ltr">
<div id="gmail-m_-5873209276696769746gmail-m_-9202678646466270310divtagdefaultwrapper" style="font-size:12pt;color:rgb(0,0,0);font-family:Calibri,Helvetica,sans-serif" dir="ltr">
<p style="margin-top:0px;margin-bottom:0px">Hello,</p>
<p style="margin-top:0px;margin-bottom:0px"><br>
</p>
<p style="margin-top:0px;margin-bottom:0px">I have a 2D system which is assembled by each elemental matrix. Ae is my elemental matrix, auxRHSe(:) and RHSe(:) and corresponding right hand side, idx is the global index. My code is as follow, however ,the assembling
rate is super slow(Marked red in the code). I am not sure whether the assembling type is right or not. Since for each element, idx are not continuous numbers. Do you have any idea what is the better way to assemble the matrix?</p>
<p style="margin-top:0px;margin-bottom:0px"><br>
</p>
<p style="margin-top:0px;margin-bottom:0px">Thanks</p>
<p style="margin-top:0px;margin-bottom:0px"></p>
<p style="margin-top:0px;margin-bottom:0px"><br>
</p>
<p style="margin-top:0px;margin-bottom:0px"></p>
<div>block</div>
<div> PetscErrorCode ierr</div>
<div> PetscMPIInt rank,nproc, mystart</div>
<div> PetscInt nelem</div>
<div> integer,allocatable ::info(:)</div>
<div> real(wp), allocatable :: Ae(:,:), auxRHSe(:),RHSe(:)</div>
<div> integer, allocatable :: idx(:)</div>
<div> PetscScalar, pointer :: xx_v(:)</div>
<div> </div>
<div> PC prec</div>
<div> PetscScalar :: val</div>
<div> Vec xvec,bvec,uvec</div>
<div> Mat Amat</div>
<div> KSP ksp</div>
<div> PetscViewer viewer</div>
<div> PetscInt geq,j,k,ne,M,Istart,Iend</div>
<div> PetscBool flg</div>
<div> KSPConvergedReason reason</div>
<div> Vec dummyVec, dummyVecs(1)</div>
<div> MatNullSpace nullspace</div>
<div><br>
</div>
<div> call PetscInitialize( PETSC_NULL_CHARACTER, ierr )</div>
<div><br>
</div>
<div> if (ierr .ne. 0) then</div>
<div> print*,'Unable to initialize PETSc'</div>
<div> stop</div>
<div> endif</div>
<div> call MPI_Comm_size(PETSC_COMM_WORLD,nproc,ierr)</div>
<div> call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)</div>
<div> mystart=rank+1</div>
<div><br>
</div>
<div> ! Parameter set</div>
<div> info=ptsystem%getInitInfo()</div>
<div> nelem = info(Info_EleMatNum)+info(Info_FixedDOFNum)+info(Info_NumConstrain) </div>
<div> print*,'nelem',nelem</div>
<div> !-------------------------------------</div>
<div> ! Create Matrix</div>
<div> call MatCreate(PETSC_COMM_WORLD,Amat,ierr)</div>
<div> call MatSetSizes( Amat,PETSC_DECIDE, PETSC_DECIDE, info(1), info(1), ierr )</div>
<div> call MatSetType( Amat, MATMPIBAIJ, ierr )</div>
<div> ! call MatMPIAIJSetPreallocation(Amat,f9,PETSC_NULL_INTEGER,f6,PETSC_NULL_INTEGER, ierr)</div>
<div> call MatSetFromOptions( Amat, ierr )</div>
<div> call MatSetUp( Amat, ierr )</div>
<div> call MatGetOwnershipRange( Amat, Istart, Iend, ierr )</div>
<div> </div>
<div> xvec = tVec(0)</div>
<div> call MatCreateVecs( Amat, PETSC_NULL_VEC, xvec, ierr )</div>
<div> call VecSetFromOptions( xvec, ierr )</div>
<div> call VecDuplicate( xvec, bvec, ierr )</div>
<div> call VecDuplicate( xvec, uvec, ierr )</div>
<div><br>
</div>
<div> t1 = MPI_WTime(); </div>
<div><br>
</div>
<div> do i=mystart,nelem,nproc</div>
<div> call ptSystem%getElementalMAT(i, Ae, auxRHSe, idx)</div>
<div> ne=size(idx)</div>
<div> if (allocated(auxRHSe)) call VecSetValues(bvec,ne,idx-1,auxRHSe,ADD_VALUES,ierr)</div>
<div> call MatSetValues(Amat,ne,idx-1,ne,idx-1,Ae,ADD_VALUES,ierr)</div>
<div> end do</div>
<div><br>
</div>
<div> nelem = info(Info_EleRHSNum)</div>
<div> mystart = rank+1</div>
<div><br>
</div>
<div> do i = mystart, nelem, nproc </div>
<div> call ptSystem%getElementalRHS(i, RHSe, idx)</div>
<div> print*,'idx',idx</div>
<div> ne=size(idx)</div>
<div> if (allocated(RHSe)) call VecSetValues(bvec,ne,idx-1,RHSe,ADD_VALUES,ierr)</div>
<div> end do </div>
<div><span style="color:rgb(255,0,0)"> call MatAssemblyBegin(Amat,MAT_FINAL_ASSEMBLY,ierr)</span></div>
<div><span style="color:rgb(255,0,0)"> call MatAssemblyEnd(Amat,MAT_FINAL_ASSEMBLY,ierr) ! this part is slow, the for loop above is done but here it may get stuck</span></div>
<div><span style="color:rgb(255,0,0)"> call VecAssemblyBegin(bvec,ierr) ! For a 2500 DOF system, assembling only takes over 2 seconds</span></div>
<div><span style="color:rgb(255,0,0)"> call VecAssemblyEnd(bvec,ierr) ! But for a 10000 DOF system , it gets stuck</span></div>
<div> t2 = MPI_WTime(); </div>
<div><span style="font-family:Calibri,Helvetica,sans-serif"> </span>print*,'assembeling time',t2-t1</div>
<div> ! Solve</div>
<div> call KSPCreate(PETSC_COMM_WORLD,ksp,ierr)</div>
<div>! Set operators. Here the matrix that defines the linear system</div>
<div>! also serves as the preconditioning matrix.</div>
<div> call KSPSetOperators(ksp,Amat,Amat,ierr)</div>
<div></div>
<div> call KSPSetFromOptions(ksp,ierr)</div>
<div>! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -</div>
<div>! Solve the linear system</div>
<div>! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -</div>
<div> call KSPSetType(ksp,KSPCG,ierr)</div>
<div> call KSPGetPC(ksp,prec,ierr)</div>
<div> ! call KSPSetPCSide(ksp,PC_SYMMETRIC,ierr)</div>
<div> call PCSetType(prec,PCJACOBI,ierr)</div>
<div> call KSPSolve(ksp,bvec,xvec,ierr)</div>
<div></div>
<div> call PetscFinalize(ierr)</div>
<div><br>
</div>
<div>end block</div>
<br>
<p></p>
<p style="margin-top:0px;margin-bottom:0px"><br>
</p>
<div id="gmail-m_-5873209276696769746gmail-m_-9202678646466270310Signature">
<div id="gmail-m_-5873209276696769746gmail-m_-9202678646466270310divtagdefaultwrapper" dir="ltr" style="font-size:12pt;color:rgb(0,0,0);background-color:rgb(255,255,255);font-family:Calibri,Arial,Helvetica,sans-serif,EmojiFont,"Apple Color Emoji","Segoe UI Emoji",NotoColorEmoji,"Segoe UI Symbol","Android Emoji",EmojiSymbols">
<p><span style="color:rgb(33,33,33);font-size:13.3333px">Yaxiong Chen, </span><br style="color:rgb(33,33,33);font-size:13.3333px">
<span style="color:rgb(33,33,33);font-size:13.3333px">Ph.D. Student </span><br style="color:rgb(33,33,33);font-size:13.3333px">
<br style="color:rgb(33,33,33);font-size:13.3333px">
<span style="font-size:13.3333px"></span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px">School of Mechanical Engineering, 3171 </span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px"><span style="color:rgb(34,34,34);font-family:Roboto,arial,sans-serif;font-size:13px">585 Purdue Mall</span></span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px">West Lafayette, IN 47907</span></p>
<div><span style="font-size:13.3333px"><br>
</span></div>
<p></p>
<p><br style="color:rgb(33,33,33);font-size:13.3333px">
<br>
</p>
</div>
</div>
</div>
</div>
</blockquote></div>