pmtst Subroutine

private subroutine pmtst(me, Ioptn, Numgr, Nparm, Param, Icntyp, Mactrk, Iact, Pttbl, Iptb, Indm, Actdif, Iphse, Iwork, Liwrk, Work, Lwrk, Confun, Pmat)

This subroutine sets up the (nparm+1) by mactrk matrix pmat. for 1 <= j <= mactrk, the top nparm elements of column j of pmat will contain the negative of the gradient of active constraint j (if constraint j is of type 2, i.e. of the form abs(f(x) - f(parwrk,x)) <= w, the left side will be treated as f(x) - f(parwrk,x) if this quantity is nonnegative and will be treated as f(parwrk,x) - f(x) otherwise). the (nparm+1)st row of pmat will contain actdif, the right side of the inequalities gradient.vector >= actdif.

Type Bound

conmax_solver

Arguments

Type IntentOptional Attributes Name
class(conmax_solver), intent(inout) :: me
integer, intent(in) :: Ioptn
integer, intent(in) :: Numgr
integer, intent(in) :: Nparm
real(kind=wp) :: Param(Nparm)
integer :: Icntyp(Numgr)
integer, intent(in) :: Mactrk
integer :: Iact(Numgr)
real(kind=wp) :: Pttbl(Iptb,Indm)
integer, intent(in) :: Iptb
integer, intent(in) :: Indm
real(kind=wp), intent(in) :: Actdif(Numgr)
integer, intent(in) :: Iphse
integer :: Iwork(Liwrk)
integer, intent(in) :: Liwrk
real(kind=wp) :: Work(Lwrk)
integer, intent(in) :: Lwrk
real(kind=wp) :: Confun(Numgr,Nparm+1)
real(kind=wp) :: Pmat(Nparm+1,Numgr)

Calls

proc~~pmtst~~CallsGraph proc~pmtst conmax_solver%pmtst proc~derst conmax_solver%derst proc~pmtst->proc~derst proc~iloc iloc proc~pmtst->proc~iloc fnset fnset proc~derst->fnset

Called by

proc~~pmtst~~CalledByGraph proc~pmtst conmax_solver%pmtst proc~rkcon conmax_solver%rkcon proc~rkcon->proc~pmtst proc~rkpar conmax_solver%rkpar proc~rkcon->proc~rkpar proc~rkpar->proc~pmtst proc~conmax conmax_solver%conmax proc~conmax->proc~rkcon

Source Code

    subroutine pmtst(me, Ioptn, Numgr, Nparm, Param, Icntyp, Mactrk, Iact, Pttbl, &
                     Iptb, Indm, Actdif, Iphse, Iwork, Liwrk, Work, Lwrk, &
                     Confun, Pmat)

        implicit none

        class(conmax_solver), intent(inout) :: me
        integer, intent(in) :: Indm
        integer, intent(in) :: Ioptn
        integer, intent(in) :: Iphse
        integer, intent(in) :: Iptb
        integer, intent(in) :: Liwrk
        integer, intent(in) :: Lwrk
        integer, intent(in) :: Mactrk
        integer, intent(in) :: Nparm
        integer, intent(in) :: Numgr
        integer :: Iwork(Liwrk)
        integer :: Iact(Numgr)
        integer :: Icntyp(Numgr)
        real(wp), intent(in) :: Actdif(Numgr)
        real(wp) :: Confun(Numgr, Nparm + 1)
        real(wp) :: Param(Nparm)
        real(wp) :: Pmat(Nparm + 1, Numgr)
        real(wp) :: Pttbl(Iptb, Indm)
        real(wp) :: Work(Lwrk)

        integer :: i, ii, ilc22, ilc24, ilc35, ioptth, ipt, j, l
        integer :: npar1

        ! SET MACHINE AND PRECISION DEPENDENT CONSTANTS FOR PMTST.
        ilc22 = iloc(22, Nparm, Numgr)
        ilc24 = iloc(24, Nparm, Numgr)
        ilc35 = iloc(35, Nparm, Numgr)
        ioptth = (Ioptn - (Ioptn/100000)*100000)/10000
        npar1 = Nparm + 1

        if (ioptth > 0) then
            ! HERE IOPTTH=1 AND WE CALL DERST TO PUT GRADIENT VALUES INTO CONFUN.
            ! IF IPHSE < 0 OR NO ICNTYP(L) IS POSITIVE, SET IPT=-1 TO TELL DERST
            ! TO COMPUTE STANDARD CONSTRAINTS ONLY, WHILE OTHERWISE SET IPT=0 TO
            ! TELL DERST TO COMPUTE ALL CONSTRAINTS.
            ipt = -1
            if (Iphse >= 0) then
                do l = 1, Numgr
                    if (Icntyp(l) > 0) then
                        ipt = 0
                        exit
                    end if
                end do
            end if
            call me%derst(Ioptn, Nparm, Numgr, Pttbl, Iptb, Indm, Param, ipt, Work(ilc24), &
                          Work(ilc35), Iwork(ilc22), Confun)
        end if

        do i = 1, Mactrk
            ii = Iact(i)
            ipt = abs(ii)

            ! HERE IOPTTH=0 AND WE HAVE NOT YET PLACED THE GRADIENT IN CONFUN, SO WE
            ! CALL DERST TO DO SO NOW.  DERST WILL ALSO COMPUTE THE
            ! CONSTRAINT VALUES, WHICH WILL NOT BE NEEDED HERE, BUT EXPECTING USERS TO
            ! WRITE FNSET SO THAT GRADIENT CALCULATIONS WILL NOT NEED FUNCTION VALUE
            ! CALCULATION RESULTS WOULD BE TOO MUCH OF A PROGRAMMING TRAP.
            if (ioptth <= 0) call me%derst(Ioptn, Nparm, Numgr, Pttbl, Iptb, Indm, &
                                           Param, ipt, Work(ilc24), Work(ilc35), &
                                           Iwork(ilc22), Confun)

            ! NOW THE GRADIENT FOR CONSTRAINT IPT IS IN CONFUN(IPT,.), AND WE PUT IT
            ! OR ITS NEGATIVE INTO PMAT.
            ! IF ICNTYP(IPT) <= 1 WE PROCEED AS IF WE HAD A -ACTIVE CONSTRAINT IN
            ! THE ICNTYP(IPT)=2 CASE.  IN ALL CASES WE PUT THE NEGATIVE OF THE
            ! CONSTRAINT GRADIENT INTO COLUMN I OF PMAT.
            if (Icntyp(ipt) > 1) then
                ! HERE ICNTYP(IPT)=2.
                if (ii > 0) then
                    ! HERE WE HAVE A +ACTIVE CONSTRAINT AT POINT IPT.
                    ! THE CONSTRAINT GRADIENT IS IN -CONFUN(IPT,.) SINCE THE LEFT SIDE OF
                    ! CONSTRAINT I IS F(X)-F(PARWRK,X) AND DERST COMPUTES THE
                    ! GRADIENT OF F(PARWRK,X).  THUS WE PUT CONFUN(IPT,.) IN COLUMN I OF PMAT.
                    do j = 1, Nparm
                        Pmat(j, i) = Confun(ipt, j + 1)
                    end do
                    cycle
                end if
            end if

            ! HERE WE HAVE A -ACTIVE TYPE 2 CONSTRAINT AT POINT -II OR AN ACTIVE
            ! CONSTRAINT OF TYPE -2, -1, OR 1 AT POINT II.
            do j = 1, Nparm
                Pmat(j, i) = -Confun(ipt, j + 1)
            end do
        end do

        ! PUT ACTDIF IN THE LAST ROW OF PMAT.
        do i = 1, Mactrk
            Pmat(npar1, i) = Actdif(i)
        end do

    end subroutine pmtst