[petsc-users] using real and complex together
Matthew Knepley
knepley at gmail.com
Mon Sep 28 19:17:28 CDT 2020
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/>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20200928/3a969fd7/attachment-0001.html>
More information about the petsc-users
mailing list