derst Subroutine

private subroutine derst(me, Ioptn, Nparm, Numgr, Pttbl, Iptb, Indm, Param, Ipt, Param1, v, Kcntyp, Confun)

This subroutine uses fnset to compute confun(i,1) and the partial derivatives of the function whose value is in confun(i,1) for certain value(s) of i. note that we do not want the icntyp computed by fnset to override the icntyp (or jcntyp) carried into this subroutine in icntyp, so we use kcntyp when we call fnset. (the icntyp computed by fnset was stored earlier through a call to ercmp1 from conmax.)

Type Bound

conmax_solver

Arguments

Type IntentOptional Attributes Name
class(conmax_solver), intent(inout) :: me
integer, intent(in) :: Ioptn
integer :: Nparm
integer :: Numgr
real(kind=wp), dimension(Iptb, Indm) :: Pttbl
integer :: Iptb
integer :: Indm
real(kind=wp), dimension(Nparm) :: Param
integer :: Ipt
real(kind=wp), dimension(Nparm) :: Param1
real(kind=wp), dimension(Numgr + 2*Nparm + 1, Nparm + 2) :: v
integer, dimension(Numgr) :: Kcntyp
real(kind=wp), dimension(Numgr, Nparm + 1) :: Confun

Calls

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

Called by

proc~~derst~~CalledByGraph proc~derst conmax_solver%derst proc~conmax conmax_solver%conmax proc~conmax->proc~derst proc~rkcon conmax_solver%rkcon proc~conmax->proc~rkcon proc~slpcon conmax_solver%slpcon proc~conmax->proc~slpcon proc~corrct conmax_solver%corrct proc~corrct->proc~derst proc~pmtst conmax_solver%pmtst proc~pmtst->proc~derst proc~rkcon->proc~derst proc~rkcon->proc~corrct proc~rkcon->proc~pmtst proc~rkpar conmax_solver%rkpar proc~rkcon->proc~rkpar proc~searsl conmax_solver%searsl proc~rkcon->proc~searsl proc~setu1 conmax_solver%setu1 proc~setu1->proc~derst proc~rkpar->proc~corrct proc~rkpar->proc~pmtst proc~searsl->proc~corrct proc~slpcon->proc~setu1 proc~slpcon->proc~searsl

Source Code

    subroutine derst(me, Ioptn, Nparm, Numgr, Pttbl, Iptb, Indm, Param, &
                     Ipt, Param1, v, Kcntyp, Confun)

        implicit none

        class(conmax_solver), intent(inout) :: me
        integer, intent(in)  :: Ioptn
        integer  :: Nparm
        integer  :: Numgr
        integer  :: Iptb
        integer  :: Indm
        integer  :: Ipt
        real(wp), dimension(Iptb, Indm)               :: Pttbl
        real(wp), dimension(Nparm)                   :: Param1
        real(wp), dimension(Nparm)                   :: Param
        real(wp), dimension(Numgr + 2*Nparm + 1, Nparm + 2) :: v
        integer, dimension(Numgr)                   :: Kcntyp
        real(wp), dimension(Numgr, Nparm + 1)           :: Confun

        real(wp) :: delt, delt2, up
        integer :: iopone, ioptth, iptkp, j, k, l, npar1

        ! IF THE ONES DIGIT OF IOPTN IS 0, WE CALL FNSET WITH INDFN=1 TO DO THE
        ! COMPUTATIONS DIRECTLY USING FORMULAS SUPPLIED BY THE USER.
        iopone = Ioptn - (Ioptn/10)*10
        if (iopone <= 0) then
            call me%fnset(Nparm, Numgr, Pttbl, Iptb, Indm, Param, Ipt, 1, Kcntyp, &
                          Confun)
            return
        else
            ! HERE THE ONES DIGIT OF IOPTN IS 1, AND WE APPROXIMATE THE PARTIAL
            ! DERIVATIVES USING CENTERED DIFFERENCE APPROXIMATIONS.
            ioptth = (Ioptn - (Ioptn/100000)*100000)/10000
            ! SET PRECISION DEPENDENT CONSTANTS.
            delt = sqrt(spcmn)
            delt2 = delt + delt
            if (ioptth <= 0) then
                ! HERE IOPONE=1 AND IOPTTH=0, AND WE WORK ONLY WITH CONSTRAINT IPT,
                ! WHERE IPT WILL BE AN INTEGER BETWEEN 1 AND NUMGR.
                ! L WILL BE THE INDEX OF THE VARIABLE WITH RESPECT TO WHICH WE ARE
                ! COMPUTING THE PARTIAL DERIVATIVE.
                do l = 1, Nparm
                    ! SET PARAM1 EQUAL TO PARAM, ECXEPT WITH ITS LTH COMPONENT INCREASED
                    ! BY DELT.
                    do j = 1, Nparm
                        Param1(j) = Param(j)
                    end do
                    Param1(l) = Param(l) + delt
                    ! NOW CALL FNSET WITH INDFN=0 TO PLACE THE FUNCTION IN CONSTRAINT
                    ! IPT EVALUATED AT POINT PARAM1 IN CONFUN(IPT,1).
                    call me%fnset(Nparm, Numgr, Pttbl, Iptb, Indm, Param1, Ipt, 0, Kcntyp, Confun)
                    up = Confun(Ipt, 1)
                    ! SET PARAM1 EQUAL TO PARAM, ECXEPT WITH ITS LTH COMOPONENT DECREASED
                    ! BY DELT, AND CALL FNSET AGAIN.
                    Param1(l) = Param(l) - delt
                    call me%fnset(Nparm, Numgr, Pttbl, Iptb, Indm, Param1, Ipt, 0, Kcntyp, Confun)
                    ! NOW WE CAN COMPUTE THE CENTERED-DIFFERENCE APPROXIMATION TO THE PARTIAL
                    ! DERIVATIVE OF THE FUNCTION IN CONSTRAINT IPT WITH RESPECT TO THE LTH
                    ! VARIABLE AT THE POINT PARAM.  THIS BELONGS IN CONFUN(IPT,L+1), AND
                    ! WE COULD PUT IT THERE NOW IF THE USER FOLLOWED DIRECTIONS AND DID NOT
                    ! CHANGE CONFUN(IPT,L+1) (SINCE INDFN=0) IN LATER FNSET CALLS, BUT TO
                    ! BE SAFE WE TEMPORARILY STORE IT IN V(L,1).
                    ! NOTE THAT V IS USED ELSEWHERE IN THE PROGRAM, BUT HERE IT IS JUST A
                    ! WORK ARRAY, WHILE THE WORK ARRAY PARAM1 IS NOT USED ELSEWHERE IN
                    ! THE PROGRAM.
                    v(l, 1) = (up - Confun(Ipt, 1))/delt2
                end do
                ! NOW COMPUTE THE VALUE OF THE FUNCTION AT PARAM, AND THEN PUT THE
                ! EARLIER-COMPUTED PARTIAL DERIVATIVES INTO CONFUN.
                call me%fnset(Nparm, Numgr, Pttbl, Iptb, Indm, Param, Ipt, 0, Kcntyp, &
                              Confun)
                do l = 1, Nparm
                    Confun(Ipt, l + 1) = v(l, 1)
                end do
                return
            else
                ! HERE IOPONE=1 AND IOPTTH=1, AND EACH TIME FNSET IS CALLED IT WILL
                ! COMPUTE VALUES FOR THE FUNCTIONS IN THE LEFT SIDES OF ALL CONSTRAINTS
                ! (EXCEPT THOSE WHERE FNSET SETS ICNTYP(I)=0) IF IPT=0, AND WILL COMPUTE
                ! VALUES FOR THE FUNCTIONS IN THE LEFT SIDES OF ALL STANDARD (I.E. TYPE
                ! -1 OR -2) CONSTRAINTS IF IPT=-1.
                ! WE FIRST SAVE IPT IN CASE THE USER CHANGES IT IN A FNSET CALL;  WE WILL
                ! RESTORE IT AFTER EACH FNSET CALL.
                iptkp = Ipt
                npar1 = Nparm + 1
                ! WE WILL COMPUTE APPROXIMATIONS TO PARTIAL DERIVATIVES FOR THOSE
                ! CONSTRAINTS WHICH FNSET IS ASKED BY IPT TO COMPUTE.  TO DETERMINE WHICH
                ! THESE ARE WE ZERO OUT KCNTYP;  AFTER A FNSET CALL, THE DESIRED
                ! CONSTRAINTS WILL BE THE CONSTRAINTS K WITH KCNTYP(K) /= 0 IF IPT=0,
                ! OR THE CONSTRAINTS K WITH KCNTYP(K) < 0 IF IPT=-1.
                do k = 1, Numgr
                    Kcntyp(k) = 0
                end do
                ! NOW FOLLOW BASICALLY THE SAME PROCEDURES AS IN THE IOPTTH=0 CASE DONE
                ! ABOVE.
                do l = 1, Nparm
                    do j = 1, Nparm
                        Param1(j) = Param(j)
                    end do
                    Param1(l) = Param(l) + delt
                    call me%fnset(Nparm, Numgr, Pttbl, Iptb, Indm, Param1, Ipt, 0, &
                                  Kcntyp, Confun)
                    Ipt = iptkp
                    do k = 1, Numgr
                        if (Ipt < 0) then
                            if (Kcntyp(k) >= 0) cycle
                        else if (Kcntyp(k) == 0) then
                            cycle
                        end if
                        ! SAVE THE UPPER NUMBERS IN COLUMN NPARM+1 OF V.
                        v(k, npar1) = Confun(k, 1)
                    end do
                    ! REVISE PARAM1 AND CALL FNSET AGAIN.
                    Param1(l) = Param(l) - delt
                    call me%fnset(Nparm, Numgr, Pttbl, Iptb, Indm, Param1, Ipt, 0, &
                                  Kcntyp, Confun)
                    Ipt = iptkp
                    do k = 1, Numgr
                        if (Ipt < 0) then
                            if (Kcntyp(k) >= 0) cycle
                        else if (Kcntyp(k) == 0) then
                            cycle
                        end if
                        ! STORE THE APPROXIMATE PARTIAL DERIVATIVES WITH RESPECT TO THE LTH
                        ! VARIABLE IN THE LTH COLUMN OF V.
                        v(k, l) = (v(k, npar1) - Confun(k, 1))/delt2
                    end do
                end do
                ! CALL FNSET AGAIN TO COMPUTE THE VALUES OF THE FUNCTIONS AT POINT
                ! PARAM, AND THEN PUT THE EARLIER-COMPUTED PARTIAL DERIVATIVES INTO
                ! CONFUN.
                call me%fnset(Nparm, Numgr, Pttbl, Iptb, Indm, Param, Ipt, 0, Kcntyp, &
                              Confun)
                do k = 1, Numgr
                    if (Ipt < 0) then
                        if (Kcntyp(k) >= 0) cycle
                    else if (Kcntyp(k) == 0) then
                        cycle
                    end if
                    do l = 1, Nparm
                        Confun(k, l + 1) = v(k, l)
                    end do
                end do
            end if
        end if

    end subroutine derst