Efficient Conversion of Local Stiffness Matrix to Global Stiffness
Hi all,
I work with FEM and was once told by a professor that assembling the stiffness matrix should be roughly 10% of computation time while the solver eats the remaining computation time. This is not the case with my code as assembly for problems larger than 10^6 DOF seems to take as long as the solver. That said, I am looking to improve the section of the code which takes values at the local stiffness level and populates the global stiffness matrix. We are using sparse storage by rows so this involves repeated searching of blocks of column information to determine where in the global array local coefficients should be stored. We use unstructured meshes so there is no order to the global node numbering but I am wondering if there is any a priori nodal renumbering?, or clever search strategy, etc? that could be employed to speed up this tedious but necessary part of the code. Cheers 
Did you allocate the pattern of the sparse matrix before usage? It sounds a bit like there are many reallocations involved. Also, you should assemble into a small dense element matrix and add the local contributions all at once to the global system to avoid as many lookups as possible. If you did both of those things, then maybe your sparse matrix implementation is in some point inefficient. Did you try different implementations, e.g. PETSc, MTL, DUNEISTL, ...?

This is what I do in my FEM code. It may be inefficient for problems
with that many DOF but it seems to work fine for up to 10^5 unknowns which is what I have tried. Also the assembly process can be made parallel fairly simply so that may be an option. ! Loop on element columns for element i do j=1,elem_dof(i) ! Global column number l=compr_map(elem_ft(i,j),1) if(l.ne.0)then ! Loop on element rows do k=1,elem_dof(i) ! Global row number m=compr_map(elem_ft(i,k),1) ! CSC format ! Loop on row pointers for global column l if(sparse_store.eq.0)then add=.false. search1:do ll=pointers(l),pointers(l+1)1 if(indices(ll).eq.m)then kk=ll add=.true. exit search1 endif enddo search1 ! CSR format ! Loop on column pointers for global column m elseif(sparse_store.eq.1)then add=.false. if(m.ne.0)then search2:do ll=irn(m),irn(m+1)1 if(jcn(ll).eq.l)then kk=ll add=.true. exit search2 endif enddo search2 endif endif if(add)then K_work(kk,1)=K_work(kk,1)+kelem(k,j) if(iassemble.eq.0.or.iassemble.eq.1)then Mass_work(kk,1)=Mass_work(kk,1)+melem(k,j) if(damp_type.eq.1)then Damp_work(kk,1)=Damp_work(kk,1)+delem(k,j) endif endif endif enddo endif enddo endif 
All times are GMT 4. The time now is 13:11. 