July 2016 - petsc-dev - mailman3.cels.anl.gov

Re: [petsc-dev] p4est configure error
by Satish Balay 07 Jul '16

07 Jul '16

On Thu, 7 Jul 2016, Tobin Isaac wrote: > On Thu, Jul 07, 2016 at 03:16:40PM +0200, Mark Adams wrote: > > p4est on Edison configures in w/o debugging, although I have to run > > configure twice. > > But with debugging it fails. I can not see any difference between the two > > configures other than debugging. > > Mark > > The log shows that you are configuring a p4est commit (3d3e1) that is > before the one I pushed specifically to avoid the fortran > configuration error (f54bb3). Please clear the git.p4est directory, Rerunning configure should use the appropriate git snapshot [assuming the appropriate self.gitcommit is listed in p4est.py]. With externalpackage git repos - 'nuking' should be rarely required.. And to use a different commit than whats listed in 'self.gitcommit' - one can use --download-PACKAGE-commit='newcommit' Or manually update/fix the repo sources and use --download-PACKAGE-commit=HEAD Satish > and confirm that you are configuring with the current 'petsc' branch > of p4est (`grep 'checking Point version' configure.log` should show > X.XXX-a9f22). > > Cheers, > Toby >

1 0

Re: [petsc-dev] p4est configure error
by Tobin Isaac 07 Jul '16

07 Jul '16

On Thu, Jul 07, 2016 at 03:16:40PM +0200, Mark Adams wrote: > p4est on Edison configures in w/o debugging, although I have to run > configure twice. > But with debugging it fails. I can not see any difference between the two > configures other than debugging. > Mark The log shows that you are configuring a p4est commit (3d3e1) that is before the one I pushed specifically to avoid the fortran configuration error (f54bb3). Please clear the git.p4est directory, and confirm that you are configuring with the current 'petsc' branch of p4est (`grep 'checking Point version' configure.log` should show X.XXX-a9f22). Cheers, Toby

1 0

p4est configure error
by Mark Adams 07 Jul '16

07 Jul '16

p4est on Edison configures in w/o debugging, although I have to run configure twice. But with debugging it fails. I can not see any difference between the two configures other than debugging. Mark

1 0

Re: [petsc-dev] Problem with assembly of Jacobian in ex12 & ex77
by Sander Arens 06 Jul '16

06 Jul '16

> Ok, that makes sense. I did not find documentation on PetscDSAddBoundary > but I managed to use the old functions to impose the correct BC. > See: http://www.mcs.anl.gov/petsc/petsc-dev/docs/manualpages/DM/PetscDSAddBounda…, but it does the same as DMAddBoundary. I must have made some changes without noticing when i tried to follow the > workings of the ex77 sourcecode, my apologies. I copied the example again > and it runs and converges fine now with the options you provided. If I > change -dim to 2 however it does not converge and / or produce the correct > Jacobian. But as far as I understand, ex77 is not intended to run in 2D. Indeed, jacobians and residuals are only implemented for 3d. I have troubles to understand what you mean with “tabulate the basis” and I > am not too familiar with the concept of dual spaces. My intention is to use > the PetscFE infrastructure to construct the global Jacobian to give to > the nonlinear solver. As a simple example, I tried using snes for solving a > linear elastic problem in 2D and PetscFE to create the Jacobian and to > impose the BC. > I do now have a program that runs and gives the correct solution if I use > the matrix-free SNES solver. -snes_type test gives a ratio of 0.45 however > and I am almost certain that the stiffness tensor I create is correct (I > checked the assembly algorithm in numpy and it gives the expected cauchy > stress tensor). The only explanation I have left is some misunderstanding > on my side about how the double contraction works in Petsc. > Ah yes, if you have a tensor like A^i_I^j_J you should switch the I and j (first index is the component of the test function, second index is the component of the trial function, third index is the derivative index of the test function and fourth index is the derivative index of the trial function). Thanks, Sander On 6 July 2016 at 15:24, <rickcha(a)googlemail.com> wrote: > > On 28 Jun 2016, at 23:45, Matthew Knepley <knepley(a)gmail.com> wrote: > > On Tue, Jun 28, 2016 at 6:47 PM, Sander Arens <Sander.Arens(a)ugent.be> > wrote: > >> Please reply to all, so that other people can also read this discussion >> on the mailing list. >> >> I do currently have problems understanding the inner workings of PetscFE. >>> For example how the element matrix is assembled. As far as I could deduce >>> from your example, g3_uu_3d is the fourth-order elasticity tensor which >>> means that we should form the double-dot product of it with the strain for >>> it to be the integrand of the gradient of the basis functions ( int_omega >>> [g3:epsilon]•grad(v) d omega). So for a linear elastic problem, according >>> to my understanding, g3 *should* be the elasticity tensor and g0,g1,g2 >>> = NULL if there is no volumetric force. And the residual should be f0 = 0 ; >>> f1= the cauchy stress tensor. >>> >>> Yes, that seems correct. >> > >> I also have troubles understanding how to apply boundary conditions. >>> There seems to be the concept of DMAddBoundary() and the >>> PetscDSSetBdJacobian() and PetscDSSetBdResidual(). I struggle to understand >>> what sets them apart from each other. My impression is that DMAddBoundary() >>> is for Dirichlet and Neumann conditions and that PetscDSSetBdJacobian() and >>> PetscDSSetBdResidual() are for BC types that do not fall in those >>> categories? >> >> Note that you can also use PetscDSAddBoundary now (if you use the master >> branch). DMAddBoundary is basically used to tell which boundary region to >> apply which bc. If it's a Dirichlet bc you have to provide a function to >> prescribe the values on the boundary, if it's a Neumann bc you just give it >> NULL. PetscDSSetBdJacobian() and PetscDSSetBdResidual() are for Neumann bcs. >> >> The reason why there's this difference for Dirichlet and Neumann bcs is >> because for finite elements you need to assemble stuff for non-zero Neumann >> bcs while for Dirichlet you set some dofs to some prescribed value. >> > Ok, that makes sense. I did not find documentation on PetscDSAddBoundary > but I managed to use the old functions to impose the correct BC. > > >> ex77 works with the options you provided me with. However the solver does >>> not seem to converge. >>> >> What do you mean that it works but it doesn't converge? Did you use all >> the options I provided you? In case not, you should send the output of >> -snes_view. >> > > The output of -ksp_view -ksp_monitor_true_residual -ksp_converged_reason > > I must have made some changes without noticing when i tried to follow the > workings of the ex77 sourcecode, my apologies. I copied the example again > and it runs and converges fine now with the options you provided. If I > change -dim to 2 however it does not converge and / or produce the correct > Jacobian. But as far as I understand, ex77 is not intended to run in 2D. > > > >> The problem I have with the documentation online so far is that I fail to >>> fully understand how to apply petscFE. I did for example not know that you >>> had to set the order of the finite element space and I have the feeling >>> that there is much more I am missing. >>> >> >> I don't really understand what you mean with "applying PetscFE". >> I think the simplest way to learn how to use it is to start from one of >> the examples and if you don't understand what a function does you look it >> up in the documentation. If you want to see what kind of options there are >> for some object, for example petscspace, you add "-help|grep petscspace" to >> the command line. >> > > PetscFE is a Ciarlet triple (P, P', K) meaning a PetscSpace for the primal > space P, PetscDualSpace for the dual space P', and reference > cell K which is a DM. The main function of PetscFE is to tabulate the > basis at quadrature points, and to act with dual space basis vectors > on functions. > > I have troubles to understand what you mean with “tabulate the basis” and > I am not too familiar with the concept of dual spaces. My intention is to > use the PetscFE infrastructure to construct the global Jacobian to give > to the nonlinear solver. As a simple example, I tried using snes for > solving a linear elastic problem in 2D and PetscFE to create the Jacobian > and to impose the BC. > I do now have a program that runs and gives the correct solution if I use > the matrix-free SNES solver. -snes_type test gives a ratio of 0.45 however > and I am almost certain that the stiffness tensor I create is correct (I > checked the assembly algorithm in numpy and it gives the expected cauchy > stress tensor). The only explanation I have left is some misunderstanding > on my side about how the double contraction works in Petsc. > > Thank you, > > M. Hartig > > Thanks, > > Matt > > >> Thanks, >> Sander >> >> On 28 June 2016 at 17:47, <rickcha(a)googlemail.com> wrote: >> >>> Yes, that was vague, I should have formulated that better. >>> What I’m trying to accomplish is writing a FEM code to calculate >>> deformation and stress of structure. As a first step analysing a linear >>> elastic and static problem in 2D, which I plan to extend to a >>> elastic-plastic transient analysis for a 3D case. >>> >>> I do currently have problems understanding the inner workings of >>> PetscFE. For example how the element matrix is assembled. As far as I could >>> deduce from your example, g3_uu_3d is the fourth-order elasticity tensor >>> which means that we should form the double-dot product of it with the >>> strain for it to be the integrand of the gradient of the basis functions ( >>> int_omega [g3:epsilon]•grad(v) d omega). So for a linear elastic problem, >>> according to my understanding, g3 *should* be the elasticity tensor and >>> g0,g1,g2 = NULL if there is no volumetric force. And the residual should be >>> f0 = 0 ; f1= the cauchy stress tensor. >>> >>> I also have troubles understanding how to apply boundary conditions. >>> There seems to be the concept of DMAddBoundary() and the >>> PetscDSSetBdJacobian() and PetscDSSetBdResidual(). I struggle to understand >>> what sets them apart from each other. My impression is that DMAddBoundary() >>> is for Dirichlet and Neumann conditions and that PetscDSSetBdJacobian() and >>> PetscDSSetBdResidual() are for BC types that do not fall in those >>> categories? >>> >>> ex77 works with the options you provided me with. However the solver >>> does not seem to converge. >>> >>> The problem I have with the documentation online so far is that I fail >>> to fully understand how to apply petscFE. I did for example not know that >>> you had to set the order of the finite element space and I have the feeling >>> that there is much more I am missing. >>> >>> Thank you for your precious time, >>> >>> M. Hartig >>> >>> >>> On 28 Jun 2016, at 14:30, Sander Arens <Sander.Arens(a)ugent.be> wrote: >>> >>> What are you trying to do with it? I don't think there's more >>> documentation on how to use PetscFE than what's already in the online >>> documentation and in the examples. >>> What do you mean with "grasp the concept"? >>> >>> Thanks, >>> Sander >>> >>> On 27 June 2016 at 17:41, <rickcha(a)googlemail.com> wrote: >>> >>>> That did do the trick with the examples, thank you. I am still >>>> struggling with my own simple example however. Is there any documentation >>>> available on how to use petscFE? I could not find anything besides the >>>> three examples. And I yet have to grasp the concept. >>>> >>>> >>>> Thanks in advance, >>>> M. Hartig >>>> >>>> On 27 Jun 2016, at 15:49, Sander Arens <Sander.Arens(a)ugent.be> wrote: >>>> >>>> I think you forgot to set the order of the finite element spaces, which >>>> default to 0. Does it work if you try for example ex77 with the following >>>> options found here? >>>> <https://bitbucket.org/petsc/petsc/src/8bc5349a152c81f25b59a9592e1e5936c0186…> >>>> >>>> Thanks, >>>> Sander >>>> >>>> On 27 June 2016 at 15:28, <rickcha(a)googlemail.com> wrote: >>>> >>>>> Dear Petsc developer team, >>>>> >>>>> while learning the ways of the petscFE I came across a curiosity I’d >>>>> like to point out to you. The Jacobians that are created from the petscFE >>>>> object seem to have only 0-entries (apart from the obviously empty regions) >>>>> in examples 12 and 77. That might be due to a limited understanding from my >>>>> side on how exactly petscFE and the mentioned examples work, but it struck >>>>> me as unusual. >>>>> >>>>> Thanks in advance, >>>>> >>>>> M. Hartig >>>> >>>> >>>> >>>> >>> >>> >> > > > -- > 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 > > >

1 0

Re: [petsc-dev] Problem with assembly of Jacobian in ex12 & ex77
by rickcha＠googlemail.com 06 Jul '16

06 Jul '16

> On 28 Jun 2016, at 23:45, Matthew Knepley <knepley(a)gmail.com> wrote: > > On Tue, Jun 28, 2016 at 6:47 PM, Sander Arens <Sander.Arens(a)ugent.be <mailto:[email protected]>> wrote: > Please reply to all, so that other people can also read this discussion on the mailing list. > > I do currently have problems understanding the inner workings of PetscFE. For example how the element matrix is assembled. As far as I could deduce from your example, g3_uu_3d is the fourth-order elasticity tensor which means that we should form the double-dot product of it with the strain for it to be the integrand of the gradient of the basis functions ( int_omega [g3:epsilon]•grad(v) d omega). So for a linear elastic problem, according to my understanding, g3 should be the elasticity tensor and g0,g1,g2 = NULL if there is no volumetric force. And the residual should be f0 = 0 ; f1= the cauchy stress tensor. > > Yes, that seems correct. > > I also have troubles understanding how to apply boundary conditions. There seems to be the concept of DMAddBoundary() and the PetscDSSetBdJacobian() and PetscDSSetBdResidual(). I struggle to understand what sets them apart from each other. My impression is that DMAddBoundary() is for Dirichlet and Neumann conditions and that PetscDSSetBdJacobian() and PetscDSSetBdResidual() are for BC types that do not fall in those categories? > Note that you can also use PetscDSAddBoundary now (if you use the master branch). DMAddBoundary is basically used to tell which boundary region to apply which bc. If it's a Dirichlet bc you have to provide a function to prescribe the values on the boundary, if it's a Neumann bc you just give it NULL. PetscDSSetBdJacobian() and PetscDSSetBdResidual() are for Neumann bcs. > The reason why there's this difference for Dirichlet and Neumann bcs is because for finite elements you need to assemble stuff for non-zero Neumann bcs while for Dirichlet you set some dofs to some prescribed value. Ok, that makes sense. I did not find documentation on PetscDSAddBoundary but I managed to use the old functions to impose the correct BC. > > ex77 works with the options you provided me with. However the solver does not seem to converge. > What do you mean that it works but it doesn't converge? Did you use all the options I provided you? In case not, you should send the output of -snes_view. > > The output of -ksp_view -ksp_monitor_true_residual -ksp_converged_reason I must have made some changes without noticing when i tried to follow the workings of the ex77 sourcecode, my apologies. I copied the example again and it runs and converges fine now with the options you provided. If I change -dim to 2 however it does not converge and / or produce the correct Jacobian. But as far as I understand, ex77 is not intended to run in 2D. > > The problem I have with the documentation online so far is that I fail to fully understand how to apply petscFE. I did for example not know that you had to set the order of the finite element space and I have the feeling that there is much more I am missing. > > I don't really understand what you mean with "applying PetscFE". > I think the simplest way to learn how to use it is to start from one of the examples and if you don't understand what a function does you look it up in the documentation. If you want to see what kind of options there are for some object, for example petscspace, you add "-help|grep petscspace" to the command line. > > PetscFE is a Ciarlet triple (P, P', K) meaning a PetscSpace for the primal space P, PetscDualSpace for the dual space P', and reference > cell K which is a DM. The main function of PetscFE is to tabulate the basis at quadrature points, and to act with dual space basis vectors > on functions. > I have troubles to understand what you mean with “tabulate the basis” and I am not too familiar with the concept of dual spaces. My intention is to use the PetscFE infrastructure to construct the global Jacobian to give to the nonlinear solver. As a simple example, I tried using snes for solving a linear elastic problem in 2D and PetscFE to create the Jacobian and to impose the BC. I do now have a program that runs and gives the correct solution if I use the matrix-free SNES solver. -snes_type test gives a ratio of 0.45 however and I am almost certain that the stiffness tensor I create is correct (I checked the assembly algorithm in numpy and it gives the expected cauchy stress tensor). The only explanation I have left is some misunderstanding on my side about how the double contraction works in Petsc. Thank you, M. Hartig > Thanks, > > Matt > > Thanks, > Sander > > On 28 June 2016 at 17:47, <rickcha(a)googlemail.com <mailto:[email protected]>> wrote: > Yes, that was vague, I should have formulated that better. > What I’m trying to accomplish is writing a FEM code to calculate deformation and stress of structure. As a first step analysing a linear elastic and static problem in 2D, which I plan to extend to a elastic-plastic transient analysis for a 3D case. > > I do currently have problems understanding the inner workings of PetscFE. For example how the element matrix is assembled. As far as I could deduce from your example, g3_uu_3d is the fourth-order elasticity tensor which means that we should form the double-dot product of it with the strain for it to be the integrand of the gradient of the basis functions ( int_omega [g3:epsilon]•grad(v) d omega). So for a linear elastic problem, according to my understanding, g3 should be the elasticity tensor and g0,g1,g2 = NULL if there is no volumetric force. And the residual should be f0 = 0 ; f1= the cauchy stress tensor. > > I also have troubles understanding how to apply boundary conditions. There seems to be the concept of DMAddBoundary() and the PetscDSSetBdJacobian() and PetscDSSetBdResidual(). I struggle to understand what sets them apart from each other. My impression is that DMAddBoundary() is for Dirichlet and Neumann conditions and that PetscDSSetBdJacobian() and PetscDSSetBdResidual() are for BC types that do not fall in those categories? > > ex77 works with the options you provided me with. However the solver does not seem to converge. > > The problem I have with the documentation online so far is that I fail to fully understand how to apply petscFE. I did for example not know that you had to set the order of the finite element space and I have the feeling that there is much more I am missing. > > Thank you for your precious time, > > M. Hartig > > >> On 28 Jun 2016, at 14:30, Sander Arens <Sander.Arens(a)ugent.be <mailto:[email protected]>> wrote: >> >> What are you trying to do with it? I don't think there's more documentation on how to use PetscFE than what's already in the online documentation and in the examples. >> What do you mean with "grasp the concept"? >> >> Thanks, >> Sander >> >> On 27 June 2016 at 17:41, <rickcha(a)googlemail.com <mailto:[email protected]>> wrote: >> That did do the trick with the examples, thank you. I am still struggling with my own simple example however. Is there any documentation available on how to use petscFE? I could not find anything besides the three examples. And I yet have to grasp the concept. >> >> >> Thanks in advance, >> M. Hartig >> >>> On 27 Jun 2016, at 15:49, Sander Arens <Sander.Arens(a)ugent.be <mailto:[email protected]>> wrote: >>> >>> I think you forgot to set the order of the finite element spaces, which default to 0. Does it work if you try for example ex77 with the following options found here? <https://bitbucket.org/petsc/petsc/src/8bc5349a152c81f25b59a9592e1e5936c0186…> >>> >>> Thanks, >>> Sander >>> >>> On 27 June 2016 at 15:28, <rickcha(a)googlemail.com <mailto:[email protected]>> wrote: >>> Dear Petsc developer team, >>> >>> while learning the ways of the petscFE I came across a curiosity I’d like to point out to you. The Jacobians that are created from the petscFE object seem to have only 0-entries (apart from the obviously empty regions) in examples 12 and 77. That might be due to a limited understanding from my side on how exactly petscFE and the mentioned examples work, but it struck me as unusual. >>> >>> Thanks in advance, >>> >>> M. Hartig >>> >> >> > > > > > > -- > 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

1 0

Get Coffee and Water for Development Titan LLC
by Jennifer Burke 05 Jul '16

05 Jul '16

Hi, My name is Jennifer at CoffeeQuotes, a service that compares quotes for Coffee, Water and Ice machines for offices. Provide your employees fresh brewed coffee, cold water and cold ice daily at the best rates available. The following information will be needed in order to compose a Quote: * Full Name * Company Name * Phone number * Zip * Coffee/Water/Ice? Finding the perfect vendor for your business won't cost a thing when filling out the above information. Best regards, Jennifer Burke Office Vendor Specialist jennifer(a)WATER4OFFICE.info www.WATER4OFFICE.info (http://www.WATER4OFFICE.info/default.aspx?c=&p=(407) 388-5337&f=&l=&z=32792&e=petsc-dev(a)mcs.anl.gov&utm_source=is&utm_medium=email&utm_campaign=ZO3) 11C Great Meadow Lane | East Hanover, NJ 07936 Click here to unsubscribe (http://www.WATER4OFFICE.info/[email protected])

1 0

Re: [petsc-dev] asm / gasm
by Mark Adams 03 Jul '16

03 Jul '16

I see, thanks, On Sun, Jul 3, 2016 at 8:00 PM, Fande Kong <fdkong.jd(a)gmail.com> wrote: > If you already have the overlapping subdomain, you do not have to set the > non-overlapping subdomain. If you already passed the overlapping > subdomain, you could not increase the overlap any more (however, we could > enable this very easily, but no many user cases, I think). > > > On Sun, Jul 3, 2016 at 10:45 AM, Mark Adams <mfadams(a)lbl.gov> wrote: > >> PetscErrorCode PCASMSetLocalSubdomains(PC pc,PetscInt n,IS is[],IS is_local[]) >>>> >>>> >>>> >>> "is" is the overlapping subdomain, and "is_local" is the subdomain >>> without overlap. >>> >>> >> So if I have some overlapping subdomains (I don't) I would have to prune >> the IS to get some non-overlapping subdomain? Is this just for something >> internal or does it affect the semantics of the solver? And I assume I can >> still set the overlap and ASM will grow my subdomains accordingly? >> > >

1 0

Re: [petsc-dev] asm / gasm
by Fande Kong 03 Jul '16

03 Jul '16

If you already have the overlapping subdomain, you do not have to set the non-overlapping subdomain. If you already passed the overlapping subdomain, you could not increase the overlap any more (however, we could enable this very easily, but no many user cases, I think). On Sun, Jul 3, 2016 at 10:45 AM, Mark Adams <mfadams(a)lbl.gov> wrote: > PetscErrorCode PCASMSetLocalSubdomains(PC pc,PetscInt n,IS is[],IS is_local[]) >>> >>> >>> >> "is" is the overlapping subdomain, and "is_local" is the subdomain >> without overlap. >> >> > So if I have some overlapping subdomains (I don't) I would have to prune > the IS to get some non-overlapping subdomain? Is this just for something > internal or does it affect the semantics of the solver? And I assume I can > still set the overlap and ASM will grow my subdomains accordingly? >

1 0

Re: [petsc-dev] asm / gasm
by Barry Smith 03 Jul '16

03 Jul '16

> On Jul 3, 2016, at 11:59 AM, Mark Adams <mfadams(a)lbl.gov> wrote: > > Thanks, I saw that. ASM converges slower than SOR :( You can try -mg_levels_pc_asm_local_type multiplicative that Matt added, but you probably need to switch to GMRES since I don't think Matt did a symmetric version of it. > > On Sun, Jul 3, 2016 at 6:50 PM, Barry Smith <bsmith(a)mcs.anl.gov> wrote: > > Here is the problem: > > if (!sz) { > IS nis; > ierr = ISCreateGeneral(PETSC_COMM_SELF, 0, NULL, PETSC_COPY_VALUES, &nis);CHKERRQ(ierr); > ierr = PCASMSetLocalSubdomains(subpc, 1, NULL, &nis);CHKERRQ(ierr); > ierr = ISDestroy(&nis);CHKERRQ(ierr); > > It was setting the first IS when it should be setting the second. When I run after this change I get reasonable results. > > 0 KSP Residual norm 2.469493285658e+02 > 1 KSP Residual norm 1.335901394170e+01 > 2 KSP Residual norm 5.623002716761e-01 > 3 KSP Residual norm 5.201650164984e-02 > 4 KSP Residual norm 3.763050190737e-03 > 5 KSP Residual norm 3.024410926322e-04 > Linear solve converged due to CONVERGED_RTOL iterations 5 > 0 KSP Residual norm 2.469493285658e-03 > 1 KSP Residual norm 1.335901394170e-04 > 2 KSP Residual norm 5.623002716761e-06 > 3 KSP Residual norm 5.201650164983e-07 > 4 KSP Residual norm 3.763050190737e-08 > 5 KSP Residual norm 3.024410926322e-09 > Linear solve converged due to CONVERGED_RTOL iterations 5 > 0 KSP Residual norm 2.469493285658e-08 > 1 KSP Residual norm 1.335901394170e-09 > 2 KSP Residual norm 5.623002716761e-11 > 3 KSP Residual norm 5.201650164984e-12 > 4 KSP Residual norm 3.763050190737e-13 > 5 KSP Residual norm 3.024410926322e-14 > > > > > On Jul 3, 2016, at 10:40 AM, Mark Adams <mfadams(a)lbl.gov> wrote: > > > > Barry, > > > > GAMG must be setting up the ASM incorrectly, on empty processors, when using the aggregates. In asm.c, I see that osm->is_local is NULL for empty processors, and this causes a problem when this branch is taken by less than 8 processes (of 8): > > > > if (osm->is_local) { // asm.c:282 > > .... > > ierr = VecScatterCreate(vec,osm->is_local[i],osm->y_local[i],isl,&osm->prolongation[i]);CHKERRQ(ierr); > > > > I get an error "MPI_Allreduce() called in different locations (code lines) on different processors ", because "vec" here is global. This gamg.c code must not be correct: > > > > if (!sz) { > > IS is; > > ierr = ISCreateGeneral(PETSC_COMM_SELF, 0, NULL, PETSC_COPY_VALUES, &is);CHKERRQ(ierr); > > ierr = PCASMSetLocalSubdomains(subpc, 1, &is, NULL);CHKERRQ(ierr); > > ierr = ISDestroy(&is);CHKERRQ(ierr); > > } else { > > .... > > } > > > > I will look into why osm->is_local seems to be NULL here when it should have an empty IS. > > > > Mark > >

1 0

Re: [petsc-dev] asm / gasm
by Mark Adams 03 Jul '16

03 Jul '16

Thanks, I saw that. ASM converges slower than SOR :( On Sun, Jul 3, 2016 at 6:50 PM, Barry Smith <bsmith(a)mcs.anl.gov> wrote: > > Here is the problem: > > if (!sz) { > IS nis; > ierr = ISCreateGeneral(PETSC_COMM_SELF, 0, NULL, > PETSC_COPY_VALUES, &nis);CHKERRQ(ierr); > ierr = PCASMSetLocalSubdomains(subpc, 1, NULL, > &nis);CHKERRQ(ierr); > ierr = ISDestroy(&nis);CHKERRQ(ierr); > > It was setting the first IS when it should be setting the second. When I > run after this change I get reasonable results. > > 0 KSP Residual norm 2.469493285658e+02 > 1 KSP Residual norm 1.335901394170e+01 > 2 KSP Residual norm 5.623002716761e-01 > 3 KSP Residual norm 5.201650164984e-02 > 4 KSP Residual norm 3.763050190737e-03 > 5 KSP Residual norm 3.024410926322e-04 > Linear solve converged due to CONVERGED_RTOL iterations 5 > 0 KSP Residual norm 2.469493285658e-03 > 1 KSP Residual norm 1.335901394170e-04 > 2 KSP Residual norm 5.623002716761e-06 > 3 KSP Residual norm 5.201650164983e-07 > 4 KSP Residual norm 3.763050190737e-08 > 5 KSP Residual norm 3.024410926322e-09 > Linear solve converged due to CONVERGED_RTOL iterations 5 > 0 KSP Residual norm 2.469493285658e-08 > 1 KSP Residual norm 1.335901394170e-09 > 2 KSP Residual norm 5.623002716761e-11 > 3 KSP Residual norm 5.201650164984e-12 > 4 KSP Residual norm 3.763050190737e-13 > 5 KSP Residual norm 3.024410926322e-14 > > > > > On Jul 3, 2016, at 10:40 AM, Mark Adams <mfadams(a)lbl.gov> wrote: > > > > Barry, > > > > GAMG must be setting up the ASM incorrectly, on empty processors, when > using the aggregates. In asm.c, I see that osm->is_local is NULL for empty > processors, and this causes a problem when this branch is taken by less > than 8 processes (of 8): > > > > if (osm->is_local) { // asm.c:282 > > .... > > ierr = > VecScatterCreate(vec,osm->is_local[i],osm->y_local[i],isl,&osm->prolongation[i]);CHKERRQ(ierr); > > > > I get an error "MPI_Allreduce() called in different locations (code > lines) on different processors ", because "vec" here is global. This gamg.c > code must not be correct: > > > > if (!sz) { > > IS is; > > ierr = ISCreateGeneral(PETSC_COMM_SELF, 0, NULL, > PETSC_COPY_VALUES, &is);CHKERRQ(ierr); > > ierr = PCASMSetLocalSubdomains(subpc, 1, &is, > NULL);CHKERRQ(ierr); > > ierr = ISDestroy(&is);CHKERRQ(ierr); > > } else { > > .... > > } > > > > I will look into why osm->is_local seems to be NULL here when it should > have an empty IS. > > > > Mark > >

1 0