On Jan 21, 2012, at 10:27 AM, Matthew Knepley wrote: > On Sat, Jan 21, 2012 at 12:19 PM, Xuefei (Rebecca) Yuan <xyuan at lbl.gov> > wrote: > Hello, > > After I was not able to get rid of the zeros in the matrix by > >>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, >>> PETSC_TRUE);CHKERRQ(ierr); > > > I would like to create the matrix by the following routine: > > > ierr = MatCreate(PETSC_COMM_WORLD, &jacobian);CHKERRQ(ierr); > ierr = MatSetType(jacobian, MATMPIAIJ);CHKERRQ(ierr); > ierr = MatSetSizes(jacobian, PETSC_DECIDE, PETSC_DECIDE, > (PetscInt)(info.mx*info.my*4), (PetscInt)(info.mx*info.my*4));CHKERRQ(ierr); > ierr = MatMPIAIJSetPreallocation(jacobian, 11, PETSC_NULL, 18, > PETSC_NULL);CHKERRQ(ierr); > > instead of using > > ierr = DMGetMatrix(DMMGGetDM(dmmg), MATAIJ, &jacobian);CHKERRQ(ierr); > > to save the memory. However, I found that the routine DMGetMatrix() was still > called inside DMMGSetSNESLocal(), therefore, the matrix I created was useless: > > What exactly are you getting out of using the DM? If it does not express the > pattern of your problem, why use it?
I do not want to use DMGetMatrix, so I use the MatCreate() to get my own matrix. But inside the call DMMGSetSNESLocal(), this DMGetMatrix was called and the matrix is set to 5776 nonzeros... > > Matt > > ********************************************* > > Breakpoint 2, main (argc=3, argv=0x7fff5fbff770) at twcartffxmhd.c:255 > 255 ierr = MatCreate(PETSC_COMM_WORLD, &jacobian);CHKERRQ(ierr); > (gdb) n > 256 ierr = MatSetType(jacobian, MATMPIAIJ);CHKERRQ(ierr); > (gdb) > 257 ierr = MatSetSizes(jacobian, PETSC_DECIDE, PETSC_DECIDE, > (PetscInt)(info.mx*info.my*4), (PetscInt)(info.mx*info.my*4));CHKERRQ(ierr); > (gdb) > 258 ierr = MatMPIAIJSetPreallocation(jacobian, 11, PETSC_NULL, 18, > PETSC_NULL);CHKERRQ(ierr); > (gdb) > > Breakpoint 3, main (argc=3, argv=0x7fff5fbff770) at twcartffxmhd.c:260 > 260 ierr = DMMGSetSNESLocal(dmmg, FormFunctionLocal, > FormJacobianLocal,0,0);CHKERRQ(ierr); > (gdb) s > DMMGSetSNESLocal_Private (dmmg=0x102004480, function=0x100016ba9 > <FormFunctionLocal>, jacobian=0x100050e83 <FormJacobianLocal>, ad_function=0, > admf_function=0) at damgsnes.c:932 > 932 PetscInt i,nlevels = dmmg[0]->nlevels; > (gdb) n > 934 PetscErrorCode > (*computejacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*) = 0; > (gdb) s > 937 PetscFunctionBegin; > (gdb) n > 938 if (jacobian) computejacobian = DMMGComputeJacobian; > (gdb) > 942 CHKMEMQ; > (gdb) > 943 ierr = PetscObjectGetCookie((PetscObject) > dmmg[0]->dm,&cookie);CHKERRQ(ierr); > (gdb) > 944 if (cookie == DM_COOKIE) { > (gdb) > 948 ierr = PetscOptionsHasName(PETSC_NULL, "-dmmg_form_function_ghost", > &flag);CHKERRQ(ierr); > (gdb) > 949 if (flag) { > (gdb) > 952 ierr = > DMMGSetSNES(dmmg,DMMGFormFunction,computejacobian);CHKERRQ(ierr); > (gdb) > Matrix Object: > type=seqaij, rows=100, cols=100 > total: nonzeros=5776, allocated nonzeros=5776 > using I-node routines: found 25 nodes, limit used is 5 > 954 for (i=0; i<nlevels; i++) { > (gdb) b twcartffxmhd.c:280 > Breakpoint 4 at 0x100004032: file twcartffxmhd.c, line 282. > (gdb) c > Continuing. > Matrix Object: > type=mpiaij, rows=100, cols=100 > total: nonzeros=0, allocated nonzeros=2900 > using I-node (on process 0) routines: found 20 nodes, limit used is 5 > > ********************************************* > > The final matrix after FormJacobianlLocal() call still have 5776 nonzeros: > > % Size = 100 100 > 2 % Nonzeros = 5776 > 3 zzz = zeros(5776,3); > > Is there a way that I can pass the matrix created via line 255-258 to > FormJacobianLocal()? > > Thanks very much! > > Cheers, > > Rebecca > > > > > > > > On Jan 20, 2012, at 5:23 PM, Matthew Knepley wrote: > >> On Fri, Jan 20, 2012 at 7:07 PM, Xuefei (Rebecca) Yuan <xyuan at lbl.gov> >> wrote: >> Here is the output: >> >> >> >> >> >> On Jan 20, 2012, at 5:01 PM, Matthew Knepley wrote: >> >>> On Fri, Jan 20, 2012 at 6:55 PM, Xuefei (Rebecca) Yuan <xyuan at lbl.gov> >>> wrote: >>> Hello Matt, >>> >>> I tried several times for 3.1-p8 and dev version by putting MatSetOption >>> >>> Are you sure your entries are exactly 0.0? >> >> >> Are you using ADD_VALUES? >> >> http://petsc.cs.iit.edu/petsc/petsc-dev/file/783e93230143/src/mat/impls/aij/seq/aij.c#l310 >> >> Matt >> >>> Matt >>> >>> 1) right after creation of the matrix: >>> >>> #ifdef petscDev >>> ierr = DMCreateMatrix(DMMGGetDM(dmmg), MATAIJ, >>> &jacobian);CHKERRQ(ierr); >>> #else >>> ierr = DAGetMatrix(DMMGGetDA(dmmg), MATAIJ, >>> &jacobian);CHKERRQ(ierr); >>> #endif >>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, >>> PETSC_TRUE);CHKERRQ(ierr); >>> >>> 2) at the beginning of the FormJacobianLocal() routine: >>> >>> PetscFunctionBegin; >>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, >>> PETSC_TRUE);CHKERRQ(ierr); >>> >>> 3) before calling MatAssemblyBegin() in FormJacobianLocal() routine: >>> >>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, >>> PETSC_TRUE);CHKERRQ(ierr); >>> ierr = MatAssemblyBegin(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>> ierr = MatAssemblyEnd(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>> >>> None of those works. What is wrong here? >>> >>> Thanks, >>> >>> R >>> >>> >>> On Jan 20, 2012, at 4:32 PM, Matthew Knepley wrote: >>> >>>> On Fri, Jan 20, 2012 at 6:28 PM, Xuefei (Rebecca) Yuan <xyuan at lbl.gov> >>>> wrote: >>>> Hello Matt, >>>> >>>> I have changed the code as >>>> >>>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, >>>> PETSC_TRUE);CHKERRQ(ierr); >>>> ierr = MatAssemblyBegin(jacobian, >>>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>> ierr = MatAssemblyEnd(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>> >>>> You have to set it before you start setting values, so we know to ignore >>>> them. >>>> >>>> Matt >>>> >>>> but still get the same result as before, the matrix still has 5776 >>>> nonzeros: >>>> >>>> % Size = 100 100 >>>> 2 % Nonzeros = 5776 >>>> 3 zzz = zeros(5776,3); >>>> >>>> Then I switch the order as >>>> >>>> ierr = MatAssemblyBegin(jacobian, >>>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>> ierr = MatAssemblyEnd(jacobian, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>> ierr = MatSetOption(jacobian, MAT_IGNORE_ZERO_ENTRIES, >>>> PETSC_TRUE);CHKERRQ(ierr); >>>> >>>> nothing changed. >>>> >>>> The version is 3.1-p8. >>>> >>>> Thanks very much! >>>> >>>> Best regards, >>>> >>>> Rebecca >>>> >>>> >>>> >>>> >>>> On Jan 20, 2012, at 4:09 PM, Matthew Knepley wrote: >>>> >>>>> On Fri, Jan 20, 2012 at 6:02 PM, Xuefei (Rebecca) Yuan <xyuan at lbl.gov> >>>>> wrote: >>>>> Hello all, >>>>> >>>>> This is a test for np=1 case of the nonzero structure of the jacobian >>>>> matrix. The jacobian matrix is created via >>>>> >>>>> ierr = >>>>> DMDACreate2d(comm,DMDA_BOUNDARY_NONE,DMDA_BOUNDARY_NONE,DMDA_STENCIL_BOX, >>>>> parameters.mxfield, parameters.myfield, PETSC_DECIDE, PETSC_DECIDE, 4, 2, >>>>> 0, 0, &da);CHKERRQ(ierr); >>>>> >>>>> ierr = DMCreateMatrix(DMMGGetDM(dmmg), MATAIJ, &jacobian);CHKERRQ(ierr); >>>>> >>>>> After creation of the jacobian matrix, >>>>> >>>>> ierr = MatAssemblyBegin(jacobian, >>>>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>>> ierr = MatAssemblyEnd(jacobian, >>>>> MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); >>>>> >>>>> PetscViewer viewer; >>>>> char fileName[120]; >>>>> sprintf(fileName, >>>>> "jacobian_after_creation.m");CHKERRQ(ierr); >>>>> >>>>> FILE * fp; >>>>> >>>>> ierr = >>>>> PetscViewerASCIIOpen(PETSC_COMM_WORLD,fileName,&viewer);CHKERRQ(ierr); >>>>> ierr = >>>>> PetscViewerSetFormat(viewer,PETSC_VIEWER_ASCII_MATLAB);CHKERRQ(ierr); >>>>> ierr = MatView (jacobian, viewer); CHKERRQ (ierr); >>>>> ierr = PetscFOpen(PETSC_COMM_WORLD,fileName,"a",&fp); >>>>> CHKERRQ(ierr); >>>>> ierr = >>>>> PetscViewerASCIIPrintf(viewer,"spy((spconvert(zzz)));\n");CHKERRQ(ierr); >>>>> ierr = PetscFClose(PETSC_COMM_WORLD,fp);CHKERRQ(ierr); >>>>> PetscViewerDestroy(&viewer); >>>>> >>>>> I took a look at the structure of the jacobian by storing it in the >>>>> matlab format, the matrix has 5776 nonzeros entries, however, those >>>>> values are all zeros at the moment as I have not insert or add any values >>>>> into it yet, the structure shows: (the following figure shows a global >>>>> replacement of 0.0 by 1.0 for those 5776 numbers) >>>>> >>>>> >>>>> >>>>> >>>>> Inside the FormJacobianLocal() function, I have selected the index to >>>>> pass to the nonzero values to jacobian, for example, >>>>> >>>>> ierr = MatSetValuesStencil(jacobian, 1, &row, 6, col, v, >>>>> INSERT_VALUES);CHKERRQ(ierr); >>>>> >>>>> where >>>>> >>>>> col[0].i = column[4].i; >>>>> col[1].i = column[5].i; >>>>> col[2].i = column[6].i; >>>>> col[3].i = column[9].i; >>>>> col[4].i = column[10].i; >>>>> col[5].i = column[12].i; >>>>> >>>>> >>>>> col[0].j = column[4].j; >>>>> col[1].j = column[5].j; >>>>> col[2].j = column[6].j; >>>>> col[3].j = column[9].j; >>>>> col[4].j = column[10].j; >>>>> col[5].j = column[12].j; >>>>> >>>>> col[0].c = column[4].c; >>>>> col[1].c = column[5].c; >>>>> col[2].c = column[6].c; >>>>> col[3].c = column[9].c; >>>>> col[4].c = column[10].c; >>>>> col[5].c = column[12].c; >>>>> >>>>> v[0] = value[4]; >>>>> v[1] = value[5]; >>>>> v[2] = value[6]; >>>>> v[3] = value[9]; >>>>> v[4] = value[10]; >>>>> v[5] = value[12]; >>>>> >>>>> and did not pass the zero entries into the jacobian matrix. However, >>>>> after inserting or adding all values to the matrix, by the same routine >>>>> above to take a look at the jacobian matrix in matlab format, the matrix >>>>> still has 5776 nonzeros, in which 1075 numbers are nonzeros, and the >>>>> other 4701 numbers are all zeros. The spy() gives >>>>> >>>>> >>>>> >>>>> >>>>> for the true nonzero structures. >>>>> >>>>> But the ksp_view will give the nonzeros number as 5776, instead of 1075: >>>>> >>>>> linear system matrix = precond matrix: >>>>> Matrix Object: Mat_0x84000000_1 1 MPI processes >>>>> type: seqaij >>>>> rows=100, cols=100 >>>>> total: nonzeros=5776, allocated nonzeros=5776 >>>>> >>>>> It is a waste of memory to have all those values of zeros been stored in >>>>> the jacobian. >>>>> >>>>> Is there anyway to get rid of those zero values in jacobian and has the >>>>> only nonzero numbers stored in jacobian? In such a case, the ksp_view >>>>> will tell that total: nonzeros=1075. >>>>> >>>>> MatSetOption(MAT_IGNORE_ZERO_ENTRIES, PETSC_TRUE); >>>>> >>>>> Matt >>>>> >>>>> Thanks very much! >>>>> >>>>> Have a nice weekend! >>>>> >>>>> Cheers, >>>>> >>>>> Rebecca >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> -- >>>>> What most experimenters take for granted before they begin their >>>>> experiments is infinitely more interesting than any results to which >>>>> their experiments lead. >>>>> -- Norbert Wiener >>>> >>>> >>>> >>>> >>>> -- >>>> What most experimenters take for granted before they begin their >>>> experiments is infinitely more interesting than any results to which their >>>> experiments lead. >>>> -- Norbert Wiener >>> >>> >>> >>> >>> -- >>> What most experimenters take for granted before they begin their >>> experiments is infinitely more interesting than any results to which their >>> experiments lead. >>> -- Norbert Wiener >> >> >> >> >> >> -- >> What most experimenters take for granted before they begin their experiments >> is infinitely more interesting than any results to which their experiments >> lead. >> -- Norbert Wiener > > > > > -- > What most experimenters take for granted before they begin their experiments > is infinitely more interesting than any results to which their experiments > lead. > -- Norbert Wiener -------------- next part -------------- An HTML attachment was scrubbed... URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20120121/29fe6855/attachment.htm>