If I load complex version of PETSc, how can Vec be real? On Monday, September 28, 2020, Matthew Knepley <[email protected]> wrote:
> On Mon, Sep 28, 2020 at 7:44 PM Sam Guo <[email protected]> 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 <[email protected]> wrote: >> >>> On Mon, Sep 28, 2020 at 6:29 PM Sam Guo <[email protected]> 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 <[email protected]> >>>> wrote: >>>> >>>>> On Mon, Sep 28, 2020 at 5:01 PM Sam Guo <[email protected]> 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 <[email protected]> >>>>>> wrote: >>>>>> >>>>>>> On Mon, Sep 28, 2020 at 3:43 PM Sam Guo <[email protected]> >>>>>>> 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 <[email protected]> >>>>>>>> wrote: >>>>>>>> >>>>>>>>> Thanks >>>>>>>>> >>>>>>>>> On Tuesday, May 26, 2020, Stefano Zampini < >>>>>>>>> [email protected]> 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 <[email protected]> >>>>>>>>>> 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 <[email protected]> wrote: >>>>>>>>>> >>>>>>>>>>> You can build PETSc with complex version, and declare some >>>>>>>>>>> variables as 'PETSC_REAL'. >>>>>>>>>>> Hong >>>>>>>>>>> >>>>>>>>>>> ------------------------------ >>>>>>>>>>> *From:* petsc-users <[email protected]> on behalf >>>>>>>>>>> of Sam Guo <[email protected]> >>>>>>>>>>> *Sent:* Tuesday, May 26, 2020 1:00 PM >>>>>>>>>>> *To:* PETSc <[email protected]> >>>>>>>>>>> *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/> >
