[petsc-users] using real and complex together
Matthew Knepley
knepley at gmail.com
Mon Sep 28 20:12:51 CDT 2020
On Mon, Sep 28, 2020 at 8:45 PM Junchao Zhang <junchao.zhang at gmail.com>
wrote:
> On Mon, Sep 28, 2020 at 7:21 PM Matthew Knepley <knepley at gmail.com> wrote:
>
>> On Mon, Sep 28, 2020 at 8:20 PM Sam Guo <sam.guo at cd-adapco.com> wrote:
>>
>>> Hi Randy and Matt,
>>> Thanks a lot. I’ll look into it.
>>>
>>
>> Another option, with less code for you but some code in PETSc, is to
>> check for a matrix with 0 imaginary part, and then
>> copy it to a real matrix and call real MUMPS. However, in this case, we
>> would have to check for real solves as well.
>>
> What do you mean? How to let mumps solve with a real matrix and a complex
> rhs?
>
Two solves.
Matt
> Thanks,
>>
>> Matt
>>
>>
>>> BR,
>>> Sam
>>>
>>> On Monday, September 28, 2020, Matthew Knepley <knepley at gmail.com>
>>> wrote:
>>>
>>>> On Mon, Sep 28, 2020 at 8:15 PM Randall Mackie <rlmackie862 at gmail.com>
>>>> wrote:
>>>>
>>>>> Sam, you can solve a complex matrix using a real version of PETSc by
>>>>> doubling the size of your matrix and spitting out real/imaginary parts.
>>>>>
>>>>> See this paper:
>>>>>
>>>>> https://epubs.siam.org/doi/abs/10.1137/S1064827500372262?mobileUi=0
>>>>>
>>>>
>>>> Thanks Randy. Yes, I meant that one.
>>>>
>>>> Matt
>>>>
>>>>
>>>>> Randy M.
>>>>>
>>>>>
>>>>> On Sep 28, 2020, at 5:12 PM, Sam Guo <sam.guo at cd-adapco.com> wrote:
>>>>>
>>>>> All I want is to solve both real and complex matrix using either real
>>>>> or complex version of PETSc. If I load complex version of PETSc, I waste
>>>>> memory solving real matrix. If I can solve complex matrix using real
>>>>> version of PETSc, I will accept it.
>>>>>
>>>>> On Monday, September 28, 2020, Sam Guo <sam.guo at cd-adapco.com> wrote:
>>>>>
>>>>>> If I load complex version of PETSc, how can Vec be real?
>>>>>>
>>>>>> On Monday, September 28, 2020, Matthew Knepley <knepley at gmail.com>
>>>>>> wrote:
>>>>>>
>>>>>>> On Mon, Sep 28, 2020 at 7:44 PM Sam Guo <sam.guo at cd-adapco.com>
>>>>>>> wrote:
>>>>>>>
>>>>>>>> I want to make sure I understand you. I load real version of
>>>>>>>> PETSc. If my input matrix is complex,
>>>>>>>
>>>>>>>
>>>>>>> You said that your matrix was real.
>>>>>>>
>>>>>>>
>>>>>>>> just get real and imagine parts of PETSc Vec
>>>>>>>
>>>>>>>
>>>>>>> No, the PETSc Vec would be real. You would have two vectors.
>>>>>>>
>>>>>>> Thanks,
>>>>>>>
>>>>>>> Matt
>>>>>>>
>>>>>>>
>>>>>>>> and do the matrix vector multiplication. Right?
>>>>>>>>
>>>>>>>> On Monday, September 28, 2020, Matthew Knepley <knepley at gmail.com>
>>>>>>>> wrote:
>>>>>>>>
>>>>>>>>> On Mon, Sep 28, 2020 at 6:29 PM Sam Guo <sam.guo at cd-adapco.com>
>>>>>>>>> wrote:
>>>>>>>>>
>>>>>>>>>> “ I think it would be much easier to just decompose your complex
>>>>>>>>>> work into real and imaginary parts and use PETSc with real scalars to
>>>>>>>>>> compute them separately.
>>>>>>>>>> Since you know your matrices have 0 imaginary part, this becomes
>>>>>>>>>> very straightforward.”
>>>>>>>>>>
>>>>>>>>>> I think this is what I am trying to do. But since I have to
>>>>>>>>>> provide matrix-vector operator(since I am using shell matrix), the Vec I
>>>>>>>>>> receive is complex. I need to convert complex vec to real one and then
>>>>>>>>>> convert it back(that’s the code I shown before).
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>> No, this is not what I am advocating. Keep your vectors real, just
>>>>>>>>> keep one for the real part and one for the imaginary part. Then you can
>>>>>>>>> just call MatMult() twice with your matrix.
>>>>>>>>>
>>>>>>>>> Thanks,
>>>>>>>>>
>>>>>>>>> Matt
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>> On Monday, September 28, 2020, Matthew Knepley <knepley at gmail.com>
>>>>>>>>>> wrote:
>>>>>>>>>>
>>>>>>>>>>> On Mon, Sep 28, 2020 at 5:01 PM Sam Guo <sam.guo at cd-adapco.com>
>>>>>>>>>>> wrote:
>>>>>>>>>>>
>>>>>>>>>>>> Hi Matt,
>>>>>>>>>>>> Since I use MUMPS as preconditioner, complex uses too much
>>>>>>>>>>>> memory if my input matrix is real. Ideally if I can compile real and
>>>>>>>>>>>> complex into different symbols (like MUMPS) , I can load both version
>>>>>>>>>>>> without conflict.
>>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> What I mean to say is that it would be great if it were as
>>>>>>>>>>> simple as using two different symbols, but unfortunately the problem is
>>>>>>>>>>> more difficult. I was trying to use
>>>>>>>>>>> the example of templates. This would be a very intrusive change
>>>>>>>>>>> no matter what technology you are using.
>>>>>>>>>>>
>>>>>>>>>>> So your main memory usage is from the MUMPS factorization, and
>>>>>>>>>>> you cannot afford to double that usage?
>>>>>>>>>>>
>>>>>>>>>>> You could consider writing a version of AIJ that stores real
>>>>>>>>>>> entries, but allows complex vector values. It would promote to complex for
>>>>>>>>>>> the row dot product.
>>>>>>>>>>> However, you would also have to do the same work for all the
>>>>>>>>>>> solves you do with MUMPS.
>>>>>>>>>>>
>>>>>>>>>>> I think it would be much easier to just decompose your complex
>>>>>>>>>>> work into real and imaginary parts and use PETSc with real scalars to
>>>>>>>>>>> compute them separately.
>>>>>>>>>>> Since you know your matrices have 0 imaginary part, this becomes
>>>>>>>>>>> very straightforward.
>>>>>>>>>>>
>>>>>>>>>>> Thanks,
>>>>>>>>>>>
>>>>>>>>>>> Matt
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>>> Thanks,
>>>>>>>>>>>> Sam
>>>>>>>>>>>>
>>>>>>>>>>>> On Mon, Sep 28, 2020 at 12:52 PM Matthew Knepley <
>>>>>>>>>>>> knepley at gmail.com> wrote:
>>>>>>>>>>>>
>>>>>>>>>>>>> On Mon, Sep 28, 2020 at 3:43 PM Sam Guo <sam.guo at cd-adapco.com>
>>>>>>>>>>>>> wrote:
>>>>>>>>>>>>>
>>>>>>>>>>>>>> Hi Stefano and PETSc dev team,
>>>>>>>>>>>>>> I want to try your suggestion to always load complex
>>>>>>>>>>>>>> version of PETSc but if my input matrix A is real, I want to create shell
>>>>>>>>>>>>>> matrix to matrix-vector and factorization using real only.
>>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> I do not think that will work as you expect. I will try to
>>>>>>>>>>>>> explain below.
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>> I still need to understand how MatRealPart works. Does it
>>>>>>>>>>>>>> just zero out the image numerical values or does it delete the image memory?
>>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> When we have complex values, we use the "complex" type to
>>>>>>>>>>>>> allocate and store them. Thus you cannot talk about just the memory to
>>>>>>>>>>>>> store imaginary parts.
>>>>>>>>>>>>> MatRealPart sets the imaginary parts of all the matrix
>>>>>>>>>>>>> elements to zero.
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>> If my input matrix A is real, how do I create a shell matrix
>>>>>>>>>>>>>> to matrix -vector multiplication y=A*x where A is real, PestcScalar =
>>>>>>>>>>>>>> complex, x and y are Vec? I notice there is a VecRealPart but it seems it
>>>>>>>>>>>>>> just zeros the image numerical values. It seems I still have to create a
>>>>>>>>>>>>>> PetscReal pointer to copy the real part of PetacScalar pointers like
>>>>>>>>>>>>>> following. Can you comment on it?
>>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> What you suggest would mean rewriting the matrix
>>>>>>>>>>>>> multiplication algorithm by hand after extracting the values. I am not sure
>>>>>>>>>>>>> if this
>>>>>>>>>>>>> is really what you want to do. Is the matrix memory really
>>>>>>>>>>>>> your limiting factor? Even if you tried to do this with templates, the
>>>>>>>>>>>>> memory
>>>>>>>>>>>>> from temporaries would be very hard to control.
>>>>>>>>>>>>>
>>>>>>>>>>>>> Thanks,
>>>>>>>>>>>>>
>>>>>>>>>>>>> Matt
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>> Thanks,
>>>>>>>>>>>>>> Sam
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> PetscScalar *px = nullptr;
>>>>>>>>>>>>>> VecGetArrayRead(x, &px);
>>>>>>>>>>>>>> PetscScalar *py = nullptr;
>>>>>>>>>>>>>> VecGetArray(y, &py);
>>>>>>>>>>>>>> int localSize = 0;
>>>>>>>>>>>>>> VecGetLocalSize(x, &localSize);
>>>>>>>>>>>>>> std::vector<PetasReal> realX(localSize); // I am using c++ to
>>>>>>>>>>>>>> call PETSc
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> //retrieve real part
>>>>>>>>>>>>>> for(int i = 0; i < localSize; i++) realX[i] =
>>>>>>>>>>>>>> PetscRealPart(px[i]);
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> // do real matrix-vector multiplication
>>>>>>>>>>>>>> // realY=A*realX
>>>>>>>>>>>>>> // here where realY is std::vector<PetscReal>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> //put real part back to py
>>>>>>>>>>>>>> for(int i = 0; i < localSize; i++) pv[i] = realY[i];
>>>>>>>>>>>>>> VecRestoreArray(y,&py);
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> On Tue, May 26, 2020 at 1:49 PM Sam Guo <
>>>>>>>>>>>>>> sam.guo at cd-adapco.com> wrote:
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> Thanks
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> On Tuesday, May 26, 2020, Stefano Zampini <
>>>>>>>>>>>>>>> stefano.zampini at gmail.com> wrote:
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> All the solvers/matrices/vectors works for PetscScalar
>>>>>>>>>>>>>>>> types (i.e. in your case complex)
>>>>>>>>>>>>>>>> If you need to solve for the real part only, you can
>>>>>>>>>>>>>>>> duplicate the matrix and call MatRealPart to zero out the imaginary part.
>>>>>>>>>>>>>>>> But the solve will always run in the complex space
>>>>>>>>>>>>>>>> You should not be worried about doubling the memory for a
>>>>>>>>>>>>>>>> matrix (i.e. real and imaginary part)
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> On May 26, 2020, at 11:28 PM, Sam Guo <
>>>>>>>>>>>>>>>> sam.guo at cd-adapco.com> wrote:
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> complex version is needed since matrix sometimes is real
>>>>>>>>>>>>>>>> and sometimes is complex. I want to solve real matrix without allocating
>>>>>>>>>>>>>>>> memory for imaginary part((except eigen pairs).
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> On Tuesday, May 26, 2020, Zhang, Hong <hzhang at mcs.anl.gov>
>>>>>>>>>>>>>>>> wrote:
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>>> You can build PETSc with complex version, and declare
>>>>>>>>>>>>>>>>> some variables as 'PETSC_REAL'.
>>>>>>>>>>>>>>>>> Hong
>>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>>> ------------------------------
>>>>>>>>>>>>>>>>> *From:* petsc-users <petsc-users-bounces at mcs.anl.gov> on
>>>>>>>>>>>>>>>>> behalf of Sam Guo <sam.guo at cd-adapco.com>
>>>>>>>>>>>>>>>>> *Sent:* Tuesday, May 26, 2020 1:00 PM
>>>>>>>>>>>>>>>>> *To:* PETSc <petsc-users at mcs.anl.gov>
>>>>>>>>>>>>>>>>> *Subject:* [petsc-users] using real and complex together
>>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>>> Dear PETSc dev team,
>>>>>>>>>>>>>>>>> Can I use both real and complex versions together?
>>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>>> Thanks,
>>>>>>>>>>>>>>>>> Sam
>>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> --
>>>>>>>>>>>>> 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://www.cse.buffalo.edu/~knepley/
>>>>>>>>>>>>> <http://www.cse.buffalo.edu/~knepley/>
>>>>>>>>>>>>>
>>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> --
>>>>>>>>>>> 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://www.cse.buffalo.edu/~knepley/
>>>>>>>>>>> <http://www.cse.buffalo.edu/~knepley/>
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>> --
>>>>>>>>> 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://www.cse.buffalo.edu/~knepley/
>>>>>>>>> <http://www.cse.buffalo.edu/~knepley/>
>>>>>>>>>
>>>>>>>>
>>>>>>>
>>>>>>> --
>>>>>>> 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://www.cse.buffalo.edu/~knepley/
>>>>>>> <http://www.cse.buffalo.edu/~knepley/>
>>>>>>>
>>>>>>
>>>>>
>>>>
>>>> --
>>>> 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://www.cse.buffalo.edu/~knepley/
>>>> <http://www.cse.buffalo.edu/~knepley/>
>>>>
>>>
>>
>> --
>> 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://www.cse.buffalo.edu/~knepley/
>> <http://www.cse.buffalo.edu/~knepley/>
>>
>
--
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://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20200928/e1dd3566/attachment.html>
More information about the petsc-users
mailing list