[petsc-users] using real and complex together
Sam Guo
sam.guo at cd-adapco.com
Mon Sep 28 18:44:34 CDT 2020
I want to make sure I understand you. I load real version of PETSc. If my
input matrix is complex, just get real and imagine parts of PETSc Vec 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/>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20200928/eb1e7dd9/attachment-0001.html>
More information about the petsc-users
mailing list