Molcas Forum

kazuma

Member

Registered: 2023-06-30

Posts: 2

The source code about preconditioned conjugate gradient (PCG) method.

Dear, all.
I know my question probably doesn't fit "Newbie Corner" but I can't find another good place.
I would like to ask you about the implementation of PCG in openmolcas.
I found the code written about the PCG in OpenMolcas/src/caspt2/pcg.f
However, three points differ from the implementation of Wikipedia (the section whose name is "The preconditioned conjugate gradient method").

https://en.wikipedia.org/wiki/Conjugate … nt_method˘

This is the code included in pcg.f

 100  CONTINUE
      CALL POVLVEC(IVECP,IVECP,OVLAPS)
      DSCALE=1.0D00/SQRT(OVLAPS(0,0))
      CALL PSCAVEC(DSCALE,IVECP,IVECP)
      CALL POVLVEC(IVECP,IVECR,OVLAPS)
      PR=OVLAPS(0,0)
      CALL SIGMA_CASPT2(1.0D00,0.0D0,IVECP,IVECT)
      CALL POVLVEC(IVECP,IVECT,OVLAPS)
      PT=OVLAPS(0,0)
      ALPHA=PR/PT
      CALL PLCVEC(ALPHA,1.0D00,IVECP,IVECX)
      CALL PLCVEC(-ALPHA,1.0D00,IVECT,IVECR)

      ...

      CALL PRESDIA(IVECR,IVECU,OVLAPS)
      UR=OVLAPS(0,0)
      BETA=PR/UR
      CALL PLCVEC(BETA,1.0D00,IVECU,IVECP)
      GOTO 100

I omit the convergence check and some other stuff (compiled in ... ).

① Why does it normalize the conjugate vector p (IVECP)?
② BETA seems the inverse of one in Wikipedia. That is,

\beta_{molcas} = 1/( \beta_{wiki} )

\beta_{molcas} = pr/ur = ( p_{k} r_{k} )/( r_{k+1}M^{-1}}r_{k+1})

\beta_{wiki} = ( r_{k+1}M^{-1}}r_{k+1}) /( r_{k}M^{-1}r_{k} )

③ The update of the conjugate vector seems different

[molcas] p_{k+1} = p_{k} + \beta_{molcas} M^{-1}r_{k+1}
[wiki]     p_{k+1} = M^{-1}r_{k+1} + \beta_{wiki} p_{k}

second and third difference doesn't seem equivalent.
If anyone knows anything, I really want to know it.

Ignacio

Administrator

From: Uppsala

Registered: 2015-11-03

Posts: 1,201

Re: The source code about preconditioned conjugate gradient (PCG) method.

Sure, p = a + b*c and p' = c + (1/b)*a are different, but proportional: p = b*p'. So after normalization both will be equivalent.

kazuma

Member

Registered: 2023-06-30

Posts: 2

Re: The source code about preconditioned conjugate gradient (PCG) method.

Wow!
I completely understand it!
I sincerely appreciate your help!