Hi, I'm new to PETSc, so that this can be a very simple question: I'm looking for something like VecGetSubVector, which I've seen it exists in the dev version but not in the released one. I need to write a smoother for a multigrid algorithm (something like a block Gauss Seidel) which can be written in matlab as for j = 1:ny P = <some matrix indices as function of j>; du(P) = L(P,P) \ ( rhs(P) - L(P,:)*du + D2(P,P)*du(P) ); end where L is a matrix (in my case the linearized Navier Stokes). I was thinking about using IS for declaring P, so that D2(P,P) can be obtained using MatGetSubMatrix. I would need the same for the vector du. Is there a way to do that without using the developer version? (I really don't feel like being "experienced with building, using and debugging PETSc). Thanks in advance Gianluca
Gianluca, The expected use is with the VecScatter object. First you create a VecScatter object with VecScatterCreate() then each time you need the "subvector" you call VecScatterBegin() followed by VecScatterEnd() Note that usually the VecScatter object is retained and used many times. Barry On Jan 3, 2011, at 5:22 AM, Gianluca Meneghello wrote:
Hi,
I'm new to PETSc, so that this can be a very simple question:
I'm looking for something like VecGetSubVector, which I've seen it exists in the dev version but not in the released one.
I need to write a smoother for a multigrid algorithm (something like a block Gauss Seidel) which can be written in matlab as
for j = 1:ny P = <some matrix indices as function of j>; du(P) = L(P,P) \ ( rhs(P) - L(P,:)*du + D2(P,P)*du(P) ); end
where L is a matrix (in my case the linearized Navier Stokes).
I was thinking about using IS for declaring P, so that D2(P,P) can be obtained using MatGetSubMatrix. I would need the same for the vector du.
Is there a way to do that without using the developer version? (I really don't feel like being "experienced with building, using and debugging PETSc).
Thanks in advance
Gianluca
Dear Barry, thanks a lot for your answer. I tried to do some experiments with MatGetSubMatrix, but I guess I'm doing something wrong as it looks like being painfully slow. I changed approach and now I'm using the ASM preconditioner. What I'm actually trying to do is to split the domain in different parts --- like interior and boundaries --- and relax (solve) each one with a different smoother (solver). In your opinion, is this the right approach? So far it looks much faster than my previous approach of extracting each submatrix. Also, please let me ask you one more thing. When using ASM with different subdomains on the same process, is the order in which the domains are solved the same as the one in which they are stored in the IS array passed to PCASMSetLocalSubdomains()? I would be interested in controlling this in order to build a downstream marching smoother. Looking at the references, I've noticed you have worked on multigrid. What I'm trying to do is close to what is described in Diskin, Thomas, Mineck, "Textbook Multigrid Efficiency for Leading Edge Stagnation", in case you already know the paper. http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/20040081104_2004082284.... Again, thanks a lot. Gianluca On 3 January 2011 17:43, Barry Smith <[email protected]> wrote:
Gianluca,
The expected use is with the VecScatter object. First you create a VecScatter object with VecScatterCreate() then each time you need the "subvector" you call VecScatterBegin() followed by VecScatterEnd() Note that usually the VecScatter object is retained and used many times.
Barry
On Jan 3, 2011, at 5:22 AM, Gianluca Meneghello wrote:
Hi,
I'm new to PETSc, so that this can be a very simple question:
I'm looking for something like VecGetSubVector, which I've seen it exists in the dev version but not in the released one.
I need to write a smoother for a multigrid algorithm (something like a block Gauss Seidel) which can be written in matlab as
for j = 1:ny P = <some matrix indices as function of j>; du(P) = L(P,P) \ ( rhs(P) - L(P,:)*du + D2(P,P)*du(P) ); end
where L is a matrix (in my case the linearized Navier Stokes).
I was thinking about using IS for declaring P, so that D2(P,P) can be obtained using MatGetSubMatrix. I would need the same for the vector du.
Is there a way to do that without using the developer version? (I really don't feel like being "experienced with building, using and debugging PETSc).
Thanks in advance
Gianluca
On Thu, Jan 6, 2011 at 04:01, Gianluca Meneghello <[email protected]>wrote:
When using ASM with different subdomains on the same process, is the order in which the domains are solved the same as the one in which they are stored in the IS array passed to PCASMSetLocalSubdomains()? I would be interested in controlling this in order to build a downstream marching smoother.
PCASM is "Additive" rather than multiplicative so there is no order. PCFieldSplit includes multiplicative and symmetric-multiplicative variants, but each subdomain ("split" in this context) is global which exposes less concurrency (subdomains can be empty on some processes so this only affects performance, not which Schwarz-like algorithms are implementable).
On Jan 6, 2011, at 6:01 AM, Gianluca Meneghello wrote:
Dear Barry,
thanks a lot for your answer.
I tried to do some experiments with MatGetSubMatrix, but I guess I'm doing something wrong as it looks like being painfully slow.
Hmm, we've always found the getsubmatrix takes a few percent of the time. Perhaps you are calling it repeatedly for the each domain, rather than once and reusing it? Also use MatGetSubMatrices() and get all the submatrices in one call rather than one at a time.
I changed approach and now I'm using the ASM preconditioner. What I'm actually trying to do is to split the domain in different parts --- like interior and boundaries --- and relax (solve) each one with a different smoother (solver). In your opinion, is this the right approach?
Worth trying since it is easy. You can experiment with different smoothers on the subdomains using the -sub_pc_type etc options and set different prefixes for different subdomains.
So far it looks much faster than my previous approach of extracting each submatrix.
ASM just uses MatGetSubMatrices() so shouldn't be faster or slower than a custom code that does the same thing.
Also, please let me ask you one more thing. When using ASM with different subdomains on the same process, is the order in which the domains are solved the same as the one in which they are stored in the IS array passed to PCASMSetLocalSubdomains()? I would be interested in controlling this in order to build a downstream marching smoother.
It is only additive, there is no order as Jed noted. Doing multiplicative in general is tricky because you want to just update the parts of the residual that need to be updated.
Looking at the references, I've noticed you have worked on multigrid. What I'm trying to do is close to what is described in Diskin, Thomas, Mineck, "Textbook Multigrid Efficiency for Leading Edge Stagnation", in case you already know the paper. http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/20040081104_2004082284....
Again, thanks a lot.
Gianluca
On 3 January 2011 17:43, Barry Smith <[email protected]> wrote:
Gianluca,
The expected use is with the VecScatter object. First you create a VecScatter object with VecScatterCreate() then each time you need the "subvector" you call VecScatterBegin() followed by VecScatterEnd() Note that usually the VecScatter object is retained and used many times.
Barry
On Jan 3, 2011, at 5:22 AM, Gianluca Meneghello wrote:
Hi,
I'm new to PETSc, so that this can be a very simple question:
I'm looking for something like VecGetSubVector, which I've seen it exists in the dev version but not in the released one.
I need to write a smoother for a multigrid algorithm (something like a block Gauss Seidel) which can be written in matlab as
for j = 1:ny P = <some matrix indices as function of j>; du(P) = L(P,P) \ ( rhs(P) - L(P,:)*du + D2(P,P)*du(P) ); end
where L is a matrix (in my case the linearized Navier Stokes).
I was thinking about using IS for declaring P, so that D2(P,P) can be obtained using MatGetSubMatrix. I would need the same for the vector du.
Is there a way to do that without using the developer version? (I really don't feel like being "experienced with building, using and debugging PETSc).
Thanks in advance
Gianluca
Dear Barry and Jed, thanks again for your answers. I'm at the moment trying to understand more about how ASM and FieldSplit works. I've started reading Barry's book "Domain Decomposition, Parallel Multilevel Methods for Elliptic Partial Differential Equations". I hope it's the right starting point. Please let me know if you have further suggested readings, taking into account that I know nothing on domain decomposition in general --- but something in multigrid. Let me ask you a couple of questions for the moment: Is PCFIELDSPLIT additive the same as PCASM (I guess no as it does not use overlap)? Or is it a Substructuring Method (I haven't yet arrived to that chapter of the book!). What's the difference between multiplicative and symmetric-multiplicative for PCFIELDSPLIT? Does symmetric-multiplicative refers to eq 1.15 of "Domain Decomposition"? Is it possible use ASM and/or FieldSplit with a matrix-free method? Thanks Gianluca Il giorno 06/gen/2011, alle ore 21.55, Barry Smith ha scritto:
On Jan 6, 2011, at 6:01 AM, Gianluca Meneghello wrote:
Dear Barry,
thanks a lot for your answer.
I tried to do some experiments with MatGetSubMatrix, but I guess I'm doing something wrong as it looks like being painfully slow.
Hmm, we've always found the getsubmatrix takes a few percent of the time. Perhaps you are calling it repeatedly for the each domain, rather than once and reusing it? Also use MatGetSubMatrices() and get all the submatrices in one call rather than one at a time.
I changed approach and now I'm using the ASM preconditioner. What I'm actually trying to do is to split the domain in different parts --- like interior and boundaries --- and relax (solve) each one with a different smoother (solver). In your opinion, is this the right approach?
Worth trying since it is easy. You can experiment with different smoothers on the subdomains using the -sub_pc_type etc options and set different prefixes for different subdomains.
So far it looks much faster than my previous approach of extracting each submatrix.
ASM just uses MatGetSubMatrices() so shouldn't be faster or slower than a custom code that does the same thing.
Also, please let me ask you one more thing. When using ASM with different subdomains on the same process, is the order in which the domains are solved the same as the one in which they are stored in the IS array passed to PCASMSetLocalSubdomains()? I would be interested in controlling this in order to build a downstream marching smoother.
It is only additive, there is no order as Jed noted. Doing multiplicative in general is tricky because you want to just update the parts of the residual that need to be updated.
Looking at the references, I've noticed you have worked on multigrid. What I'm trying to do is close to what is described in Diskin, Thomas, Mineck, "Textbook Multigrid Efficiency for Leading Edge Stagnation", in case you already know the paper. http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/20040081104_2004082284....
Again, thanks a lot.
Gianluca
On 3 January 2011 17:43, Barry Smith <[email protected]> wrote:
Gianluca,
The expected use is with the VecScatter object. First you create a VecScatter object with VecScatterCreate() then each time you need the "subvector" you call VecScatterBegin() followed by VecScatterEnd() Note that usually the VecScatter object is retained and used many times.
Barry
On Jan 3, 2011, at 5:22 AM, Gianluca Meneghello wrote:
Hi,
I'm new to PETSc, so that this can be a very simple question:
I'm looking for something like VecGetSubVector, which I've seen it exists in the dev version but not in the released one.
I need to write a smoother for a multigrid algorithm (something like a block Gauss Seidel) which can be written in matlab as
for j = 1:ny P = <some matrix indices as function of j>; du(P) = L(P,P) \ ( rhs(P) - L(P,:)*du + D2(P,P)*du(P) ); end
where L is a matrix (in my case the linearized Navier Stokes).
I was thinking about using IS for declaring P, so that D2(P,P) can be obtained using MatGetSubMatrix. I would need the same for the vector du.
Is there a way to do that without using the developer version? (I really don't feel like being "experienced with building, using and debugging PETSc).
Thanks in advance
Gianluca
On Jan 13, 2011, at 12:00 PM, Gianluca Meneghello wrote:
Dear Barry and Jed,
thanks again for your answers. I'm at the moment trying to understand more about how ASM and FieldSplit works. I've started reading Barry's book "Domain Decomposition, Parallel Multilevel Methods for Elliptic Partial Differential Equations". I hope it's the right starting point.
Please let me know if you have further suggested readings, taking into account that I know nothing on domain decomposition in general --- but something in multigrid.
Let me ask you a couple of questions for the moment: Is PCFIELDSPLIT additive the same as PCASM (I guess no as it does not use overlap)?
It can be additive or multiplicative depending on what you set with PCFieldSplitSetType()
Or is it a Substructuring Method (I haven't yet arrived to that chapter of the book!).
No, nothing to do with substructuring.
What's the difference between multiplicative and symmetric-multiplicative for PCFIELDSPLIT? Does symmetric-multiplicative refers to eq 1.15 of "Domain Decomposition"?
Yes
Is it possible use ASM and/or FieldSplit with a matrix-free method?
Yes, BUT the algorithms are coded around MatGetSubMatrix() and or MatGetSubMatrices() so to do matrix free you need to have code that applies "part" of the operator at a time (that is you cannot just have a matrix vector product that applies the entire operator to the entire vector. Once you have the ability to apply "part" of the operator at a time you need to code up a MATSHELL that responds appropriately to MatGetSubMatrix() and or MatGetSubMatrices() and returns new matrix-free shell matrices that apply only "their" part of the operator. This is non-trivial for many people but possible. Barry
Thanks
Gianluca
Il giorno 06/gen/2011, alle ore 21.55, Barry Smith ha scritto:
On Jan 6, 2011, at 6:01 AM, Gianluca Meneghello wrote:
Dear Barry,
thanks a lot for your answer.
I tried to do some experiments with MatGetSubMatrix, but I guess I'm doing something wrong as it looks like being painfully slow.
Hmm, we've always found the getsubmatrix takes a few percent of the time. Perhaps you are calling it repeatedly for the each domain, rather than once and reusing it? Also use MatGetSubMatrices() and get all the submatrices in one call rather than one at a time.
I changed approach and now I'm using the ASM preconditioner. What I'm actually trying to do is to split the domain in different parts --- like interior and boundaries --- and relax (solve) each one with a different smoother (solver). In your opinion, is this the right approach?
Worth trying since it is easy. You can experiment with different smoothers on the subdomains using the -sub_pc_type etc options and set different prefixes for different subdomains.
So far it looks much faster than my previous approach of extracting each submatrix.
ASM just uses MatGetSubMatrices() so shouldn't be faster or slower than a custom code that does the same thing.
Also, please let me ask you one more thing. When using ASM with different subdomains on the same process, is the order in which the domains are solved the same as the one in which they are stored in the IS array passed to PCASMSetLocalSubdomains()? I would be interested in controlling this in order to build a downstream marching smoother.
It is only additive, there is no order as Jed noted. Doing multiplicative in general is tricky because you want to just update the parts of the residual that need to be updated.
Looking at the references, I've noticed you have worked on multigrid. What I'm trying to do is close to what is described in Diskin, Thomas, Mineck, "Textbook Multigrid Efficiency for Leading Edge Stagnation", in case you already know the paper. http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/20040081104_2004082284....
Again, thanks a lot.
Gianluca
On 3 January 2011 17:43, Barry Smith <[email protected]> wrote:
Gianluca,
The expected use is with the VecScatter object. First you create a VecScatter object with VecScatterCreate() then each time you need the "subvector" you call VecScatterBegin() followed by VecScatterEnd() Note that usually the VecScatter object is retained and used many times.
Barry
On Jan 3, 2011, at 5:22 AM, Gianluca Meneghello wrote:
Hi,
I'm new to PETSc, so that this can be a very simple question:
I'm looking for something like VecGetSubVector, which I've seen it exists in the dev version but not in the released one.
I need to write a smoother for a multigrid algorithm (something like a block Gauss Seidel) which can be written in matlab as
for j = 1:ny P = <some matrix indices as function of j>; du(P) = L(P,P) \ ( rhs(P) - L(P,:)*du + D2(P,P)*du(P) ); end
where L is a matrix (in my case the linearized Navier Stokes).
I was thinking about using IS for declaring P, so that D2(P,P) can be obtained using MatGetSubMatrix. I would need the same for the vector du.
Is there a way to do that without using the developer version? (I really don't feel like being "experienced with building, using and debugging PETSc).
Thanks in advance
Gianluca
Barry, I'm afraid I didn't explain myself well in my first question. If I use PCFIELDSPLIT with PCFieldSplitSetType(pc,PC_COMPOSITE_ADDITIVE), is it the same as using PCASM? Concerning the MATSHELL to be used with PCFIELDSPLIT, is there an example from where to start from? I guess I'm one of the people for which it is non-trivial :-) Thanks Gianluca Il giorno 13/gen/2011, alle ore 20.14, Barry Smith ha scritto:
On Jan 13, 2011, at 12:00 PM, Gianluca Meneghello wrote:
Dear Barry and Jed,
thanks again for your answers. I'm at the moment trying to understand more about how ASM and FieldSplit works. I've started reading Barry's book "Domain Decomposition, Parallel Multilevel Methods for Elliptic Partial Differential Equations". I hope it's the right starting point.
Please let me know if you have further suggested readings, taking into account that I know nothing on domain decomposition in general --- but something in multigrid.
Let me ask you a couple of questions for the moment: Is PCFIELDSPLIT additive the same as PCASM (I guess no as it does not use overlap)?
It can be additive or multiplicative depending on what you set with PCFieldSplitSetType()
Or is it a Substructuring Method (I haven't yet arrived to that chapter of the book!).
No, nothing to do with substructuring.
What's the difference between multiplicative and symmetric-multiplicative for PCFIELDSPLIT? Does symmetric-multiplicative refers to eq 1.15 of "Domain Decomposition"?
Yes
Is it possible use ASM and/or FieldSplit with a matrix-free method?
Yes, BUT the algorithms are coded around MatGetSubMatrix() and or MatGetSubMatrices() so to do matrix free you need to have code that applies "part" of the operator at a time (that is you cannot just have a matrix vector product that applies the entire operator to the entire vector. Once you have the ability to apply "part" of the operator at a time you need to code up a MATSHELL that responds appropriately to MatGetSubMatrix() and or MatGetSubMatrices() and returns new matrix-free shell matrices that apply only "their" part of the operator. This is non-trivial for many people but possible.
Barry
Thanks
Gianluca
Il giorno 06/gen/2011, alle ore 21.55, Barry Smith ha scritto:
On Jan 6, 2011, at 6:01 AM, Gianluca Meneghello wrote:
Dear Barry,
thanks a lot for your answer.
I tried to do some experiments with MatGetSubMatrix, but I guess I'm doing something wrong as it looks like being painfully slow.
Hmm, we've always found the getsubmatrix takes a few percent of the time. Perhaps you are calling it repeatedly for the each domain, rather than once and reusing it? Also use MatGetSubMatrices() and get all the submatrices in one call rather than one at a time.
I changed approach and now I'm using the ASM preconditioner. What I'm actually trying to do is to split the domain in different parts --- like interior and boundaries --- and relax (solve) each one with a different smoother (solver). In your opinion, is this the right approach?
Worth trying since it is easy. You can experiment with different smoothers on the subdomains using the -sub_pc_type etc options and set different prefixes for different subdomains.
So far it looks much faster than my previous approach of extracting each submatrix.
ASM just uses MatGetSubMatrices() so shouldn't be faster or slower than a custom code that does the same thing.
Also, please let me ask you one more thing. When using ASM with different subdomains on the same process, is the order in which the domains are solved the same as the one in which they are stored in the IS array passed to PCASMSetLocalSubdomains()? I would be interested in controlling this in order to build a downstream marching smoother.
It is only additive, there is no order as Jed noted. Doing multiplicative in general is tricky because you want to just update the parts of the residual that need to be updated.
Looking at the references, I've noticed you have worked on multigrid. What I'm trying to do is close to what is described in Diskin, Thomas, Mineck, "Textbook Multigrid Efficiency for Leading Edge Stagnation", in case you already know the paper. http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/20040081104_2004082284....
Again, thanks a lot.
Gianluca
On 3 January 2011 17:43, Barry Smith <[email protected]> wrote:
Gianluca,
The expected use is with the VecScatter object. First you create a VecScatter object with VecScatterCreate() then each time you need the "subvector" you call VecScatterBegin() followed by VecScatterEnd() Note that usually the VecScatter object is retained and used many times.
Barry
On Jan 3, 2011, at 5:22 AM, Gianluca Meneghello wrote:
Hi,
I'm new to PETSc, so that this can be a very simple question:
I'm looking for something like VecGetSubVector, which I've seen it exists in the dev version but not in the released one.
I need to write a smoother for a multigrid algorithm (something like a block Gauss Seidel) which can be written in matlab as
for j = 1:ny P = <some matrix indices as function of j>; du(P) = L(P,P) \ ( rhs(P) - L(P,:)*du + D2(P,P)*du(P) ); end
where L is a matrix (in my case the linearized Navier Stokes).
I was thinking about using IS for declaring P, so that D2(P,P) can be obtained using MatGetSubMatrix. I would need the same for the vector du.
Is there a way to do that without using the developer version? (I really don't feel like being "experienced with building, using and debugging PETSc).
Thanks in advance
Gianluca
On Jan 13, 2011, at 2:09 PM, Gianluca Meneghello wrote:
Barry,
I'm afraid I didn't explain myself well in my first question.
If I use PCFIELDSPLIT with PCFieldSplitSetType(pc,PC_COMPOSITE_ADDITIVE), is it the same as using PCASM?
Same in what sense? It solves a bunch of subproblems independently and adds together all the solutions. There can be overlapping in the fields or not depending how what you choose. The decomposition in ASM is by "geometry" while the decomposition in the PCFIELDSPLIT is between different "fields" or "types of variables". So yes they have many similarities.
Concerning the MATSHELL to be used with PCFIELDSPLIT, is there an example from where to start from? I guess I'm one of the people for which it is non-trivial :-)
Not really. Barry
Thanks
Gianluca
Il giorno 13/gen/2011, alle ore 20.14, Barry Smith ha scritto:
On Jan 13, 2011, at 12:00 PM, Gianluca Meneghello wrote:
Dear Barry and Jed,
thanks again for your answers. I'm at the moment trying to understand more about how ASM and FieldSplit works. I've started reading Barry's book "Domain Decomposition, Parallel Multilevel Methods for Elliptic Partial Differential Equations". I hope it's the right starting point.
Please let me know if you have further suggested readings, taking into account that I know nothing on domain decomposition in general --- but something in multigrid.
Let me ask you a couple of questions for the moment: Is PCFIELDSPLIT additive the same as PCASM (I guess no as it does not use overlap)?
It can be additive or multiplicative depending on what you set with PCFieldSplitSetType()
Or is it a Substructuring Method (I haven't yet arrived to that chapter of the book!).
No, nothing to do with substructuring.
What's the difference between multiplicative and symmetric-multiplicative for PCFIELDSPLIT? Does symmetric-multiplicative refers to eq 1.15 of "Domain Decomposition"?
Yes
Is it possible use ASM and/or FieldSplit with a matrix-free method?
Yes, BUT the algorithms are coded around MatGetSubMatrix() and or MatGetSubMatrices() so to do matrix free you need to have code that applies "part" of the operator at a time (that is you cannot just have a matrix vector product that applies the entire operator to the entire vector. Once you have the ability to apply "part" of the operator at a time you need to code up a MATSHELL that responds appropriately to MatGetSubMatrix() and or MatGetSubMatrices() and returns new matrix-free shell matrices that apply only "their" part of the operator. This is non-trivial for many people but possible.
Barry
Thanks
Gianluca
Il giorno 06/gen/2011, alle ore 21.55, Barry Smith ha scritto:
On Jan 6, 2011, at 6:01 AM, Gianluca Meneghello wrote:
Dear Barry,
thanks a lot for your answer.
I tried to do some experiments with MatGetSubMatrix, but I guess I'm doing something wrong as it looks like being painfully slow.
Hmm, we've always found the getsubmatrix takes a few percent of the time. Perhaps you are calling it repeatedly for the each domain, rather than once and reusing it? Also use MatGetSubMatrices() and get all the submatrices in one call rather than one at a time.
I changed approach and now I'm using the ASM preconditioner. What I'm actually trying to do is to split the domain in different parts --- like interior and boundaries --- and relax (solve) each one with a different smoother (solver). In your opinion, is this the right approach?
Worth trying since it is easy. You can experiment with different smoothers on the subdomains using the -sub_pc_type etc options and set different prefixes for different subdomains.
So far it looks much faster than my previous approach of extracting each submatrix.
ASM just uses MatGetSubMatrices() so shouldn't be faster or slower than a custom code that does the same thing.
Also, please let me ask you one more thing. When using ASM with different subdomains on the same process, is the order in which the domains are solved the same as the one in which they are stored in the IS array passed to PCASMSetLocalSubdomains()? I would be interested in controlling this in order to build a downstream marching smoother.
It is only additive, there is no order as Jed noted. Doing multiplicative in general is tricky because you want to just update the parts of the residual that need to be updated.
Looking at the references, I've noticed you have worked on multigrid. What I'm trying to do is close to what is described in Diskin, Thomas, Mineck, "Textbook Multigrid Efficiency for Leading Edge Stagnation", in case you already know the paper. http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/20040081104_2004082284....
Again, thanks a lot.
Gianluca
On 3 January 2011 17:43, Barry Smith <[email protected]> wrote:
Gianluca,
The expected use is with the VecScatter object. First you create a VecScatter object with VecScatterCreate() then each time you need the "subvector" you call VecScatterBegin() followed by VecScatterEnd() Note that usually the VecScatter object is retained and used many times.
Barry
On Jan 3, 2011, at 5:22 AM, Gianluca Meneghello wrote:
Hi,
I'm new to PETSc, so that this can be a very simple question:
I'm looking for something like VecGetSubVector, which I've seen it exists in the dev version but not in the released one.
I need to write a smoother for a multigrid algorithm (something like a block Gauss Seidel) which can be written in matlab as
for j = 1:ny P = <some matrix indices as function of j>; du(P) = L(P,P) \ ( rhs(P) - L(P,:)*du + D2(P,P)*du(P) ); end
where L is a matrix (in my case the linearized Navier Stokes).
I was thinking about using IS for declaring P, so that D2(P,P) can be obtained using MatGetSubMatrix. I would need the same for the vector du.
Is there a way to do that without using the developer version? (I really don't feel like being "experienced with building, using and debugging PETSc).
Thanks in advance
Gianluca
Dear Barry In order to use matrix-free methods with PCFIELDSPLIT, is it correct to start from MatShellSetOperation(mat,MATOP_MATGETSUBMATRIX, (void(*)(void)) PetscErrorCode (*UserMatGetSubMatrix)(......)); or does it require to dive deeper inside PETSc code? Thanks Gianluca On 13 January 2011 22:16, Barry Smith <[email protected]> wrote:
On Jan 13, 2011, at 2:09 PM, Gianluca Meneghello wrote:
Barry,
I'm afraid I didn't explain myself well in my first question.
If I use PCFIELDSPLIT with PCFieldSplitSetType(pc,PC_COMPOSITE_ADDITIVE), is it the same as using PCASM?
Same in what sense? It solves a bunch of subproblems independently and adds together all the solutions. There can be overlapping in the fields or not depending how what you choose. The decomposition in ASM is by "geometry" while the decomposition in the PCFIELDSPLIT is between different "fields" or "types of variables". So yes they have many similarities.
Concerning the MATSHELL to be used with PCFIELDSPLIT, is there an example from where to start from? I guess I'm one of the people for which it is non-trivial :-)
Not really.
Barry
Thanks
Gianluca
Il giorno 13/gen/2011, alle ore 20.14, Barry Smith ha scritto:
On Jan 13, 2011, at 12:00 PM, Gianluca Meneghello wrote:
Dear Barry and Jed,
thanks again for your answers. I'm at the moment trying to understand more about how ASM and FieldSplit works. I've started reading Barry's book "Domain Decomposition, Parallel Multilevel Methods for Elliptic Partial Differential Equations". I hope it's the right starting point.
Please let me know if you have further suggested readings, taking into account that I know nothing on domain decomposition in general --- but something in multigrid.
Let me ask you a couple of questions for the moment: Is PCFIELDSPLIT additive the same as PCASM (I guess no as it does not use overlap)?
It can be additive or multiplicative depending on what you set with PCFieldSplitSetType()
Or is it a Substructuring Method (I haven't yet arrived to that chapter of the book!).
No, nothing to do with substructuring.
What's the difference between multiplicative and symmetric-multiplicative for PCFIELDSPLIT? Does symmetric-multiplicative refers to eq 1.15 of "Domain Decomposition"?
Yes
Is it possible use ASM and/or FieldSplit with a matrix-free method?
Yes, BUT the algorithms are coded around MatGetSubMatrix() and or MatGetSubMatrices() so to do matrix free you need to have code that applies "part" of the operator at a time (that is you cannot just have a matrix vector product that applies the entire operator to the entire vector. Once you have the ability to apply "part" of the operator at a time you need to code up a MATSHELL that responds appropriately to MatGetSubMatrix() and or MatGetSubMatrices() and returns new matrix-free shell matrices that apply only "their" part of the operator. This is non-trivial for many people but possible.
Barry
Thanks
Gianluca
Il giorno 06/gen/2011, alle ore 21.55, Barry Smith ha scritto:
On Jan 6, 2011, at 6:01 AM, Gianluca Meneghello wrote:
Dear Barry,
thanks a lot for your answer.
I tried to do some experiments with MatGetSubMatrix, but I guess I'm doing something wrong as it looks like being painfully slow.
Hmm, we've always found the getsubmatrix takes a few percent of the time. Perhaps you are calling it repeatedly for the each domain, rather than once and reusing it? Also use MatGetSubMatrices() and get all the submatrices in one call rather than one at a time.
I changed approach and now I'm using the ASM preconditioner. What I'm actually trying to do is to split the domain in different parts --- like interior and boundaries --- and relax (solve) each one with a different smoother (solver). In your opinion, is this the right approach?
Worth trying since it is easy. You can experiment with different smoothers on the subdomains using the -sub_pc_type etc options and set different prefixes for different subdomains.
So far it looks much faster than my previous approach of extracting each submatrix.
ASM just uses MatGetSubMatrices() so shouldn't be faster or slower than a custom code that does the same thing.
Also, please let me ask you one more thing. When using ASM with different subdomains on the same process, is the order in which the domains are solved the same as the one in which they are stored in the IS array passed to PCASMSetLocalSubdomains()? I would be interested in controlling this in order to build a downstream marching smoother.
It is only additive, there is no order as Jed noted. Doing multiplicative in general is tricky because you want to just update the parts of the residual that need to be updated.
Looking at the references, I've noticed you have worked on multigrid. What I'm trying to do is close to what is described in Diskin, Thomas, Mineck, "Textbook Multigrid Efficiency for Leading Edge Stagnation", in case you already know the paper. http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/20040081104_2004082284....
Again, thanks a lot.
Gianluca
On 3 January 2011 17:43, Barry Smith <[email protected]> wrote:
Gianluca,
The expected use is with the VecScatter object. First you create a VecScatter object with VecScatterCreate() then each time you need the "subvector" you call VecScatterBegin() followed by VecScatterEnd() Note that usually the VecScatter object is retained and used many times.
Barry
On Jan 3, 2011, at 5:22 AM, Gianluca Meneghello wrote:
Hi,
I'm new to PETSc, so that this can be a very simple question:
I'm looking for something like VecGetSubVector, which I've seen it exists in the dev version but not in the released one.
I need to write a smoother for a multigrid algorithm (something like a block Gauss Seidel) which can be written in matlab as
for j = 1:ny P = <some matrix indices as function of j>; du(P) = L(P,P) \ ( rhs(P) - L(P,:)*du + D2(P,P)*du(P) ); end
where L is a matrix (in my case the linearized Navier Stokes).
I was thinking about using IS for declaring P, so that D2(P,P) can be obtained using MatGetSubMatrix. I would need the same for the vector du.
Is there a way to do that without using the developer version? (I really don't feel like being "experienced with building, using and debugging PETSc).
Thanks in advance
Gianluca
On Sat, Jan 15, 2011 at 8:23 AM, Gianluca Meneghello <[email protected]>wrote:
Dear Barry
In order to use matrix-free methods with PCFIELDSPLIT, is it correct to start from
MatShellSetOperation(mat,MATOP_MATGETSUBMATRIX, (void(*)(void)) PetscErrorCode (*UserMatGetSubMatrix)(......));
or does it require to dive deeper inside PETSc code?
That is fine. You must return a Mat object, probably a MatShell in your case, which applys the requested block of the operator. Matt
Thanks
Gianluca
On 13 January 2011 22:16, Barry Smith <[email protected]> wrote:
On Jan 13, 2011, at 2:09 PM, Gianluca Meneghello wrote:
Barry,
I'm afraid I didn't explain myself well in my first question.
If I use PCFIELDSPLIT with
PCFieldSplitSetType(pc,PC_COMPOSITE_ADDITIVE), is it the same as using PCASM?
Same in what sense? It solves a bunch of subproblems independently and
adds together all the solutions. There can be overlapping in the fields or not depending how what you choose. The decomposition in ASM is by "geometry" while the decomposition in the PCFIELDSPLIT is between different "fields" or "types of variables". So yes they have many similarities.
Concerning the MATSHELL to be used with PCFIELDSPLIT, is there an
example from where to start from? I guess I'm one of the people for which it is non-trivial :-)
Not really.
Barry
Thanks
Gianluca
Il giorno 13/gen/2011, alle ore 20.14, Barry Smith ha scritto:
On Jan 13, 2011, at 12:00 PM, Gianluca Meneghello wrote:
Dear Barry and Jed,
thanks again for your answers. I'm at the moment trying to understand
more about how ASM and FieldSplit works. I've started reading Barry's book "Domain Decomposition, Parallel Multilevel Methods for Elliptic Partial Differential Equations". I hope it's the right starting point.
Please let me know if you have further suggested readings, taking into
account that I know nothing on domain decomposition in general --- but something in multigrid.
Let me ask you a couple of questions for the moment: Is PCFIELDSPLIT additive the same as PCASM (I guess no as it does not
use overlap)?
It can be additive or multiplicative depending on what you set with PCFieldSplitSetType()
Or is it a Substructuring Method (I haven't yet arrived to that chapter of the book!).
No, nothing to do with substructuring.
What's the difference between multiplicative and symmetric-multiplicative for PCFIELDSPLIT? Does symmetric-multiplicative refers to eq 1.15 of "Domain Decomposition"?
Yes
Is it possible use ASM and/or FieldSplit with a matrix-free method?
Yes, BUT the algorithms are coded around MatGetSubMatrix() and or MatGetSubMatrices() so to do matrix free you need to have code that applies "part" of the operator at a time (that is you cannot just have a matrix vector product that applies the entire operator to the entire vector. Once you have the ability to apply "part" of the operator at a time you need to code up a MATSHELL that responds appropriately to MatGetSubMatrix() and or MatGetSubMatrices() and returns new matrix-free shell matrices that apply only "their" part of the operator. This is non-trivial for many people but possible.
Barry
Thanks
Gianluca
Il giorno 06/gen/2011, alle ore 21.55, Barry Smith ha scritto:
On Jan 6, 2011, at 6:01 AM, Gianluca Meneghello wrote:
Dear Barry,
thanks a lot for your answer.
I tried to do some experiments with MatGetSubMatrix, but I guess I'm doing something wrong as it looks like being painfully slow.
Hmm, we've always found the getsubmatrix takes a few percent of the
time. Perhaps you are calling it repeatedly for the each domain, rather than once and reusing it? Also use MatGetSubMatrices() and get all the submatrices in one call rather than one at a time.
I changed approach and now I'm using the ASM preconditioner. What
I'm
actually trying to do is to split the domain in different parts --- like interior and boundaries --- and relax (solve) each one with a different smoother (solver). In your opinion, is this the right approach?
Worth trying since it is easy. You can experiment with different smoothers on the subdomains using the -sub_pc_type etc options and set different prefixes for different subdomains.
So far it looks much faster than my previous approach of extracting each submatrix.
ASM just uses MatGetSubMatrices() so shouldn't be faster or slower than a custom code that does the same thing.
Also, please let me ask you one more thing. When using ASM with different subdomains on the same process, is the order in which the domains are solved the same as the one in which they are stored in
the
IS array passed to PCASMSetLocalSubdomains()? I would be interested in controlling this in order to build a downstream marching smoother.
It is only additive, there is no order as Jed noted. Doing multiplicative in general is tricky because you want to just update the parts of the residual that need to be updated.
Looking at the references, I've noticed you have worked on
multigrid.
What I'm trying to do is close to what is described in Diskin, Thomas, Mineck, "Textbook Multigrid Efficiency for Leading Edge Stagnation", in case you already know the paper.
http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/20040081104_2004082284....
Again, thanks a lot.
Gianluca
On 3 January 2011 17:43, Barry Smith <[email protected]> wrote:
Gianluca,
The expected use is with the VecScatter object. First you create a
VecScatter object with VecScatterCreate() then each time you need the "subvector" you call VecScatterBegin() followed by VecScatterEnd() Note that usually the VecScatter object is retained and used many times.
Barry
On Jan 3, 2011, at 5:22 AM, Gianluca Meneghello wrote:
Hi,
I'm new to PETSc, so that this can be a very simple question:
I'm looking for something like VecGetSubVector, which I've seen it exists in the dev version but not in the released one.
I need to write a smoother for a multigrid algorithm (something
like a
block Gauss Seidel) which can be written in matlab as
for j = 1:ny P = <some matrix indices as function of j>; du(P) = L(P,P) \ ( rhs(P) - L(P,:)*du + D2(P,P)*du(P) ); end
where L is a matrix (in my case the linearized Navier Stokes).
I was thinking about using IS for declaring P, so that D2(P,P) can be obtained using MatGetSubMatrix. I would need the same for the vector du.
Is there a way to do that without using the developer version? (I really don't feel like being "experienced with building, using and debugging PETSc).
Thanks in advance
Gianluca
-- 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
On Thu, Jan 13, 2011 at 18:16, Barry Smith <[email protected]> wrote:
If I use PCFIELDSPLIT with PCFieldSplitSetType(pc,PC_COMPOSITE_ADDITIVE), is it the same as using PCASM?
Same in what sense? It solves a bunch of subproblems independently and adds together all the solutions. There can be overlapping in the fields or not depending how what you choose. The decomposition in ASM is by "geometry" while the decomposition in the PCFIELDSPLIT is between different "fields" or "types of variables". So yes they have many similarities.
To put it differently, PCASM exposes as much concurrency as possible, making it efficient to use with spatially "local" subdomains. In contrast, PCFIELDSPLIT decomposes a multi-physics problem into sub-problems that are hopefully "better understood" so that efficient solvers are available. While you can make PCFIELDSPLIT use the same decomposition as PCASM uses, so that the algorithm is functionally equivalent, it does not expose the same concurrency, so would not scale in parallel.
participants (4)
-
Barry Smith -
Gianluca Meneghello -
Jed Brown -
Matthew Knepley