Ok! thanks! It runs so good!! I supposed it when I read yours fisrts mails but I tried it last day and I didn't write the code correctly. When I have readed this mail I've tried again and I've found the error!
jordi En/na Matthew Knepley ha escrit: > The code you have there will create a matrix for every loop iteration, > and will only destroy one at the end. The new MatMatMult() interface > has an argument which allows the user to pass in a matrix. > > Matt > > On 5/10/06, *Jordi Marc? Nogu?* <jordi.marce at upc.edu > <mailto:jordi.marce at upc.edu>> wrote: > > Yes, of course I destroy the matrix.... The scheme of my code when I > create a diagonal Mass-lumping matrix is the code below. > > In MatSeqAIJSetPreallocation(M_aux1,6,PETSC_NULL); the number is 6 > because in a general coordinates it's possible this matrix changes in a > 6x6 full matrix. > > "element3D *p_element = lfiber[i].getElement(e)" is a internal > procedure to obtain information about the element > > ------------------------------------------------------------------ > > Mat M_aux1, M_aux2; > > MatCreateSeqAIJ(PETSC_COMM_SELF, 6, 6, 6, PETSC_NULL, &M_aux1); > MatCreateSeqAIJ(PETSC_COMM_SELF, 6, 6, 6, PETSC_NULL, &M_aux2); > > MatSeqAIJSetPreallocation(M_aux1,6,PETSC_NULL); > MatSeqAIJSetPreallocation(M_aux2,6,PETSC_NULL); > > MatSetFromOptions(M_aux1); > MatSetOption(M_aux1, MAT_SYMMETRIC); > MatSetOption(M_aux1, MAT_IGNORE_ZERO_ENTRIES); > > MatSetFromOptions(M_aux2); > MatSetOption(M_aux2, MAT_SYMMETRIC); > MatSetOption(M_aux2, MAT_IGNORE_ZERO_ENTRIES); > > > for(uint32_t i=0; i<nfiber;i++) > { > for(uint32_t e=0; e<lfiber[i].nelements;e++) > { > > MatZeroentries(M_aux1); > MatZeroentries(M_aux2); > > element3D *p_element = lfiber[i].getElement(e); > > p_element->updateTs(); > > MatSetValue(M_aux1, 0, 0, p_element->L0 * value.pho / 6., > INSERT_VALUES); > MatSetValue(M_aux1, 1, 1, p_element->L0 * value.pho / 6., > INSERT_VALUES); > MatSetValue(M_aux1, 2, 2, p_element->L0 * value.pho / 6., > INSERT_VALUES); > MatSetValue(M_aux1, 3, 3, p_element->L0 * value.pho / 6., > INSERT_VALUES); > MatSetValue(M_aux1, 4, 4, p_element->L0 * value.pho / 6., > INSERT_VALUES); > MatSetValue(M_aux1, 5, 5, p_element->L0 * value.pho / 6., > INSERT_VALUES); > > > MatAssemblyBegin(M_aux1,MAT_FINAL_ASSEMBLY); > MatAssemblyEnd(M_aux1,MAT_FINAL_ASSEMBLY); > > > MatMatMult(M_aux1,p_element->T,&M_aux2); > MatMatMult(p_element->TT,M_aux2,&M_aux1); > > } > > // here I work with the matrix M_aux2, but in this point the > memory > // is constant. The code doesn't waste memory > } > > MatDestroy(M_aux1); > MatDestroy(M_aux2); > > ------------------------ > > > Thanks, > > best regards > > > -- > Jordi Marc?-Nogu? > Dept. Resist?ncia de Materials i Estructures a l'Enginyeria > Universitat Polit?cnica de Catalunya (UPC) > > Edifici T45 - despatx 137 > ETSEIAT (Terrassa) > > phone: +34 937 398 728 > mail: jordi.marce at upc.edu <mailto:jordi.marce at upc.edu> > > > > > -- > "Failure has a thousand explanations. Success doesn't need one" -- Sir > Alec Guiness -- Jordi Marc?-Nogu? Dept. Resist?ncia de Materials i Estructures a l'Enginyeria Universitat Polit?cnica de Catalunya (UPC) Edifici T45 - despatx 137 ETSEIAT (Terrassa) phone: +34 937 398 728 mail: jordi.marce at upc.edu
