[petsc-users] using real and complex together
Sam Guo
sam.guo at cd-adapco.com
Mon Sep 28 19:12:47 CDT 2020
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/>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20200928/7efffe77/attachment-0001.html>
More information about the petsc-users
mailing list