On Fri, Jan 24, 2025 at 10:36 AM MIGUEL MOLINOS PEREZ <mmoli...@us.es> wrote:
> Thanks Matt, I tried that too, and the problem remains. The box is updated > only if I set no periodic bcc. > What do you mean by "The box is updated"? I am trying to understand how you test things. Clearly the coordinates are updated, even in the periodic case. Thus, I need to understand the test. Once we do that, we can work backwards to the first broken thing. Thanks, Matt > Miguel > > On 24 Jan 2025, at 14:20, Matthew Knepley <knep...@gmail.com> wrote: > > On Fri, Jan 24, 2025 at 4:41 AM MIGUEL MOLINOS PEREZ <mmoli...@us.es> > wrote: > >> Dear Matt, the error was in the implementation of the volume expansion >> function. I updated it, and it works finte under finite domains. However, >> if I include periodic boundary conditions the volume of the cell does not >> accommodate the volume expansion of the particles. The deformation gradient >> is not the identity… I guess I am missing the fine detail on how periodic >> bcc are implemented in DMDA mesh, I’m right? >> > > DMDA identifies vertices using a VecScatter to implement periodic BC. This > should be insensitive to coordinates. However, I don't think the algorithm > below is correct for local coordinates. You use GlobalToLocal(), which > means that some global coordinate "wins" for each local cell, so cells on > the periodic boundary can be wrong. I would set the local coordinates by > hand as well. > > Thanks, > > Matt > > >> Thanks, >> Miguel >> >> static PetscErrorCode Volumetric_Expansion_DMDA(DM * da, >> const Eigen::Matrix3d& F) { >> >> PetscInt i, j, mstart, m, nstart, n, pstart, p, k; >> Vec local, global; >> DMDACoor3d ***coors, ***coorslocal; >> DM cda; >> >> PetscFunctionBeginUser; >> PetscCall(DMGetCoordinateDM(*da, &cda)); >> PetscCall(DMGetCoordinates(*da, &global)); >> PetscCall(DMGetCoordinatesLocal(*da, &local)); >> PetscCall(DMDAVecGetArray(cda, global, &coors)); >> PetscCall(DMDAVecGetArrayRead(cda, local, &coorslocal)); >> PetscCall(DMDAGetCorners(cda, &mstart, &nstart, &pstart, &m, &n, &p)); >> for (i = mstart; i < mstart + m; i++) { >> for (j = nstart; j < nstart + n; j++) { >> for (k = pstart; k < pstart + p; k++) { >> coors[k][j][i].x = coorslocal[k][j][i].x * F(0, 0); >> coors[k][j][i].y = coorslocal[k][j][i].y * F(1, 1); >> coors[k][j][i].z = coorslocal[k][j][i].z * F(2, 2); >> } >> } >> } >> PetscCall(DMDAVecRestoreArray(cda, global, &coors)); >> PetscCall(DMDAVecRestoreArrayRead(cda, local, &coorslocal)); >> >> PetscCall(DMGlobalToLocalBegin(cda, global, INSERT_VALUES, local)); >> PetscCall(DMGlobalToLocalEnd(cda, global, INSERT_VALUES, local)); >> >> PetscFunctionReturn(PETSC_SUCCESS); >> } >> >> On 17 Jan 2025, at 18:01, MIGUEL MOLINOS PEREZ <mmoli...@us.es> wrote: >> >> You are right!! Thank you again! >> >> Miguel >> >> On Jan 17, 2025, at 5:18 PM, Matthew Knepley <knep...@gmail.com> wrote: >> >> On Fri, Jan 17, 2025 at 10:49 AM MIGUEL MOLINOS PEREZ <mmoli...@us.es> >> wrote: >> >>> Now the error is in the call to DMSwarmMigrate >>> >> >> You have almost certainly overwritten memory somewhere. Can you use >> vlagrind or Address Sanitizer? >> >> Thanks, >> >> Matt >> >> >>> Miguel >>> >>> [0]PETSC ERROR: >>> ------------------------------------------------------------------------ >>> [0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, >>> probably memory access out of range >>> [0]PETSC ERROR: Try option -start_in_debugger or >>> -on_error_attach_debugger >>> [0]PETSC ERROR: or see >>> https://urldefense.us/v3/__https://petsc.org/release/faq/*valgrind__;Iw!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSupuIaVdEg$ >>> and >>> https://urldefense.us/v3/__https://petsc.org/release/faq/__;!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSupsTWdVgG$ >>> >>> [0]PETSC ERROR: --------------------- Stack Frames >>> ------------------------------------ >>> [0]PETSC ERROR: The line numbers in the error traceback are not always >>> exact. >>> [0]PETSC ERROR: #1 DMSwarmDataBucketGetSizes() at >>> /Users/migmolper/petsc/src/dm/impls/swarm/data_bucket.c:297 >>> [0]PETSC ERROR: #2 DMSwarmMigrate_CellDMScatter() at >>> /Users/migmolper/petsc/src/dm/impls/swarm/swarm_migrate.c:201 >>> [0]PETSC ERROR: #3 DMSwarmMigrate() at >>> /Users/migmolper/petsc/src/dm/impls/swarm/swarm.c:1349 >>> [0]PETSC ERROR: #4 main() at >>> /Users/migmolper/DMD/driver-tasting-SOLERA.cpp:41 >>> application called MPI_Abort(MPI_COMM_WORLD, 59) - process 0 >>> >>> On Jan 17, 2025, at 4:22 PM, Matthew Knepley <knep...@gmail.com> wrote: >>> >>> On Fri, Jan 17, 2025 at 10:08 AM MIGUEL MOLINOS PEREZ <mmoli...@us.es> >>> wrote: >>> >>>> Thank you Matt, this the piece of code I use to change the coordinates >>>> of the DM obtained using: >>>> >>> >>> You do not need the call to DMSetCoordinates(). What happens when you >>> remove it? >>> >>> Thanks, >>> >>> Matt >>> >>> >>>> >>>> DMSwarmGetCellDM(Simulation.atomistic_data, &bounding_cell); >>>> DMGetApplicationContext(bounding_cell, &background_mesh); >>>> >>>> Thanks, >>>> Miguel >>>> >>>> >>>> /************************************************************************/ >>>> >>>> PetscErrorCode Volumetric_Expansion(DM dm, const Eigen::Matrix3d& F) { >>>> PetscErrorCode ierr; >>>> Vec coordinates; >>>> PetscScalar* coordArray; >>>> PetscInt xs, ys, zs, xm, ym, zm, i, j, k; >>>> PetscInt dim, M, N, P; >>>> >>>> PetscFunctionBegin; >>>> // Get DMDA information >>>> ierr = DMDAGetInfo(dm, &dim, &M, &N, &P, NULL, NULL, NULL, NULL, NULL, >>>> NULL, >>>> NULL, NULL, NULL); >>>> CHKERRQ(ierr); >>>> ierr = DMDAGetCorners(dm, &xs, &ys, &zs, &xm, &ym, &zm); >>>> CHKERRQ(ierr); >>>> >>>> // Get the coordinates vector >>>> ierr = DMGetCoordinates(dm, &coordinates); >>>> CHKERRQ(ierr); >>>> ierr = VecGetArray(coordinates, &coordArray); >>>> CHKERRQ(ierr); >>>> >>>> // Update the coordinates based on the desired transformation >>>> for (k = zs; k < zs + zm; k++) { >>>> for (j = ys; j < ys + ym; j++) { >>>> for (i = xs; i < xs + xm; i++) { >>>> PetscInt idx = >>>> ((k * N + j) * M + i) * dim; // Index for the i, j, k point >>>> coordArray[idx] = coordArray[idx] * F(0,0); // Update x-coordinate >>>> coordArray[idx + 1] = coordArray[idx + 1] * F(1,1); // Update >>>> y-coordinate >>>> coordArray[idx + 2] = coordArray[idx + 2] * F(2,2); // Update >>>> z-coordinate >>>> } >>>> } >>>> } >>>> >>>> // Restore the coordinates vector >>>> ierr = VecRestoreArray(coordinates, &coordArray); >>>> CHKERRQ(ierr); >>>> >>>> // Set the updated coordinates back to the DMDA >>>> ierr = DMSetCoordinates(dm, coordinates); >>>> CHKERRQ(ierr); >>>> >>>> PetscFunctionReturn(0); >>>> } >>>> >>>> >>>> /************************************************************************/ >>>> >>>> On 17 Jan 2025, at 16:00, Matthew Knepley <knep...@gmail.com> wrote: >>>> >>>> On Fri, Jan 17, 2025 at 9:45 AM MIGUEL MOLINOS PEREZ <mmoli...@us.es> >>>> wrote: >>>> >>>>> I tried what you suggested, but still I got this error message. Maybe >>>>> I should use main release? >>>>> >>>> >>>> No. I suspect something is wrong with the way you are setting >>>> coordinates. Can you share the code? >>>> >>>> Thanks, >>>> >>>> Matt >>>> >>>> >>>>> Miguel >>>>> >>>>> [4]PETSC ERROR: >>>>> ------------------------------------------------------------------------ >>>>> [4]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, >>>>> probably memory access out of range >>>>> [4]PETSC ERROR: Try option -start_in_debugger or >>>>> -on_error_attach_debugger >>>>> [4]PETSC ERROR: or see >>>>> https://urldefense.us/v3/__https://petsc.org/release/faq/*valgrind__;Iw!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSupuIaVdEg$ >>>>> and >>>>> https://urldefense.us/v3/__https://petsc.org/release/faq/__;!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSupsTWdVgG$ >>>>> >>>>> [4]PETSC ERROR: --------------------- Stack Frames >>>>> ------------------------------------ >>>>> [4]PETSC ERROR: The line numbers in the error traceback are not always >>>>> exact. >>>>> [4]PETSC ERROR: #1 Pack_PetscReal_1_0() at >>>>> /Users/migmolper/petsc/src/vec/is/sf/impls/basic/sfpack.c:373 >>>>> [4]PETSC ERROR: #2 PetscSFLinkPackRootData_Private() at >>>>> /Users/migmolper/petsc/src/vec/is/sf/impls/basic/sfpack.c:932 >>>>> [4]PETSC ERROR: #3 PetscSFLinkPackRootData() at >>>>> /Users/migmolper/petsc/src/vec/is/sf/impls/basic/sfpack.c:966 >>>>> [4]PETSC ERROR: #4 PetscSFBcastBegin_Basic() at >>>>> /Users/migmolper/petsc/src/vec/is/sf/impls/basic/sfbasic.c:357 >>>>> [4]PETSC ERROR: #5 PetscSFBcastWithMemTypeBegin() at >>>>> /Users/migmolper/petsc/src/vec/is/sf/interface/sf.c:1513 >>>>> [4]PETSC ERROR: #6 VecScatterBegin_Internal() at >>>>> /Users/migmolper/petsc/src/vec/is/sf/interface/vscat.c:70 >>>>> [4]PETSC ERROR: #7 VecScatterBegin() at >>>>> /Users/migmolper/petsc/src/vec/is/sf/interface/vscat.c:1316 >>>>> [4]PETSC ERROR: #8 DMGlobalToLocalBegin_DA() at >>>>> /Users/migmolper/petsc/src/dm/impls/da/dagtol.c:15 >>>>> [4]PETSC ERROR: #9 DMGlobalToLocalBegin() at >>>>> /Users/migmolper/petsc/src/dm/interface/dm.c:2844 >>>>> [4]PETSC ERROR: #10 DMGetCoordinatesLocalSetUp() at >>>>> /Users/migmolper/petsc/src/dm/interface/dmcoordinates.c:565 >>>>> [4]PETSC ERROR: #11 DMGetCoordinatesLocal() at >>>>> /Users/migmolper/petsc/src/dm/interface/dmcoordinates.c:599 >>>>> [4]PETSC ERROR: #12 _DMLocatePoints_DMDARegular_IS() at >>>>> /Users/migmolper/DMD/SOLERA/Atoms/Atom.cpp:531 >>>>> [4]PETSC ERROR: #13 DMLocatePoints_DMDARegular() at >>>>> /Users/migmolper/DMD/SOLERA/Atoms/Atom.cpp:586 >>>>> [4]PETSC ERROR: #14 DMLocatePoints() at >>>>> /Users/migmolper/petsc/src/dm/interface/dmcoordinates.c:1194 >>>>> [4]PETSC ERROR: #15 DMSwarmMigrate_CellDMScatter() at >>>>> /Users/migmolper/petsc/src/dm/impls/swarm/swarm_migrate.c:219 >>>>> [4]PETSC ERROR: #16 DMSwarmMigrate() at >>>>> /Users/migmolper/petsc/src/dm/impls/swarm/swarm.c:1349 >>>>> [4]PETSC ERROR: #17 main() at >>>>> /Users/migmolper/DMD/driver-tasting-SOLERA.cpp:41 >>>>> >>>>> >>>>> >>>>> On Jan 15, 2025, at 4:56 PM, MIGUEL MOLINOS PEREZ <mmoli...@us.es> >>>>> wrote: >>>>> >>>>> Thank you Matt for the useful info. I’ll try your idea. >>>>> >>>>> Miguel >>>>> >>>>> On 15 Jan 2025, at 16:48, Matthew Knepley <knep...@gmail.com> wrote: >>>>> >>>>> On Wed, Jan 15, 2025 at 10:41 AM MIGUEL MOLINOS PEREZ <mmoli...@us.es> >>>>> wrote: >>>>> >>>>>> Thank you Matt. >>>>>> >>>>>> Yes, I am getting the "CellDM" from the DMSwarm. >>>>>> >>>>>> 1. I have recently overhauled this functionality because it was not >>>>>> flexible enough for the plasma simulation we do. Thus main and release >>>>>> work >>>>>> differently. >>>>>> >>>>>> >>>>>> Nice to hear that. Should I move to main? >>>>>> >>>>> >>>>> The changes allow you to have several cell DMs. I want to bin >>>>> particles in space, but also in velocity, and then in the tensor product >>>>> of >>>>> space and velocity. Moreover, sometimes I want to use different Swarm >>>>> fields as the DM field for the solver. You can do all that with main now. >>>>> If you just need a single DM with the same DM fields, release is fine. >>>>> >>>>> >>>>>> 2. I assume you are using release >>>>>> >>>>>> >>>>>> You are correct. >>>>>> >>>>>> 3. In both main and release, if you change the coordinates of your >>>>>> CellDM mesh, you need to rebin the particles. The easiest way to do this >>>>>> is >>>>>> to call DMSwarmMigrate(sw, PETSC_FALSE). >>>>>> >>>>>> >>>>>> What do you mean by rebin? >>>>>> >>>>> >>>>> When you provide the cell DM, Swrm makes a "sort context" that bins >>>>> the particles into DM cells. If you change the coordinates, this binning >>>>> will change, so you need it to "rebin" or recreate the sort context. >>>>> >>>>> Thanks, >>>>> >>>>> Matt >>>>> >>>>> >>>>>> Miguel >>>>>> >>>>>> >>>>>> Thanks, >>>>>> >>>>>> Matt >>>>>> >>>>>> >>>>>>> Best, >>>>>>> Miguel >>>>>>> >>>>>> -- >>>>>> 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 >>>>>> >>>>>> https://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSuphgzDu4D$ >>>>>> >>>>>> <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSupiN9pla9$ >>>>>> > >>>>>> >>>>>> >>>>>> >>>>> >>>>> -- >>>>> 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 >>>>> >>>>> https://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSuphgzDu4D$ >>>>> >>>>> <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSupiN9pla9$ >>>>> > >>>>> >>>>> >>>>> >>>>> >>>> >>>> -- >>>> 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 >>>> >>>> https://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSuphgzDu4D$ >>>> >>>> <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSupiN9pla9$ >>>> > >>>> >>>> >>>> >>> >>> -- >>> 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 >>> >>> https://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSuphgzDu4D$ >>> >>> <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSupiN9pla9$ >>> > >>> >>> >>> >> >> -- >> 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 >> >> https://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSuphgzDu4D$ >> >> <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSupiN9pla9$ >> > >> >> >> >> > > -- > 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 > > https://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSuphgzDu4D$ > > <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSupiN9pla9$ > > > > > -- 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 https://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSuphgzDu4D$ <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!bl1IKWYD5vUTh2-bmOUYT72lASrjcNt_e-FTAgWDrFKIB9bal_DAR9Nx9MMl8oqM4OEMlkElerSupiN9pla9$ >