example4 Program

Uses

  • program~~example4~~UsesGraph program~example4 example4 module~example4_model example4_model program~example4->module~example4_model module~odrpack odrpack program~example4->module~odrpack module~odrpack_kinds odrpack_kinds program~example4->module~odrpack_kinds module~example4_model->module~odrpack_kinds module~odrpack->module~odrpack_kinds iso_fortran_env iso_fortran_env module~odrpack->iso_fortran_env module~odrpack_kinds->iso_fortran_env

Default ODR job, with parameter bounds. This sample problem comes from Zwolak et al. 2001 (High Performance Computing Symposium, "Estimating rate constants in cell cycle models"). The call to ODRPACK95 is modified from the call the authors make to ODRPACK. This is done to illustrate the need for bounds. The authors could just have easily used the call statement here to solve their problem. Curious users are encouraged to remove the bounds in the call statement, run the code, and compare the results to the current call statement.


Calls

program~~example4~~CallsGraph program~example4 example4 proc~odr odr program~example4->proc~odr proc~dodcnt dodcnt proc~odr->proc~dodcnt proc~dodpe1 dodpe1 proc~odr->proc~dodpe1 proc~dodphd dodphd proc~odr->proc~dodphd proc~workspace_dimensions workspace_dimensions proc~odr->proc~workspace_dimensions proc~doddrv doddrv proc~dodcnt->proc~doddrv interface~dcopy dcopy proc~doddrv->interface~dcopy interface~ddot ddot proc~doddrv->interface~ddot interface~dnrm2 dnrm2 proc~doddrv->interface~dnrm2 proc~derstep derstep proc~doddrv->proc~derstep proc~detaf detaf proc~doddrv->proc~detaf proc~dfctrw dfctrw proc~doddrv->proc~dfctrw proc~dflags dflags proc~doddrv->proc~dflags proc~diniwk diniwk proc~doddrv->proc~diniwk proc~diwinf diwinf proc~doddrv->proc~diwinf proc~djck djck proc~doddrv->proc~djck proc~dodchk dodchk proc~doddrv->proc~dodchk proc~dodmn dodmn proc~doddrv->proc~dodmn proc~dodper dodper proc~doddrv->proc~dodper proc~dpack dpack proc~doddrv->proc~dpack proc~dsetn dsetn proc~doddrv->proc~dsetn proc~dunpac dunpac proc~doddrv->proc~dunpac proc~dwght dwght proc~doddrv->proc~dwght proc~dwinf dwinf proc~doddrv->proc~dwinf proc~mbfb mbfb proc~doddrv->proc~mbfb proc~dhstep dhstep proc~derstep->proc~dhstep proc~dfctr dfctr proc~dfctrw->proc~dfctr proc~diniwk->interface~dcopy proc~diniwk->proc~dflags proc~dsclb dsclb proc~diniwk->proc~dsclb proc~dscld dscld proc~diniwk->proc~dscld proc~djck->proc~dhstep proc~djckm djckm proc~djck->proc~djckm proc~dodmn->interface~dcopy proc~dodmn->interface~ddot proc~dodmn->interface~dnrm2 proc~dodmn->proc~dflags proc~dodmn->proc~dpack proc~dodmn->proc~dunpac proc~dodmn->proc~dwght proc~dacces dacces proc~dodmn->proc~dacces proc~devjac devjac proc~dodmn->proc~devjac proc~dodlm dodlm proc~dodmn->proc~dodlm proc~dodpcr dodpcr proc~dodmn->proc~dodpcr proc~dodvcv dodvcv proc~dodmn->proc~dodvcv proc~dodper->proc~dodpe1 proc~dodper->proc~dodphd proc~dodpe2 dodpe2 proc~dodper->proc~dodpe2 proc~dodpe3 dodpe3 proc~dodper->proc~dodpe3 proc~dpack->interface~dcopy proc~dunpac->interface~dcopy proc~mbfb->proc~dhstep proc~dacces->proc~diwinf proc~dacces->proc~dwinf proc~devjac->interface~ddot proc~devjac->proc~dunpac proc~devjac->proc~dwght proc~difix difix proc~devjac->proc~difix proc~djaccd djaccd proc~devjac->proc~djaccd proc~djacfd djacfd proc~devjac->proc~djacfd proc~dfctr->interface~ddot proc~djckc djckc proc~djckm->proc~djckc proc~djckz djckz proc~djckm->proc~djckz proc~dpvb dpvb proc~djckm->proc~dpvb proc~dpvd dpvd proc~djckm->proc~dpvd proc~dodlm->interface~ddot proc~dodlm->interface~dnrm2 proc~dodlm->proc~dwght proc~dodstp dodstp proc~dodlm->proc~dodstp proc~dscale dscale proc~dodlm->proc~dscale proc~dodpcr->proc~dodphd proc~dodpcr->proc~dflags proc~dodpc1 dodpc1 proc~dodpcr->proc~dodpc1 proc~dodpc2 dodpc2 proc~dodpcr->proc~dodpc2 proc~dodpc3 dodpc3 proc~dodpcr->proc~dodpc3 dpodi dpodi proc~dodvcv->dpodi proc~dodvcv->proc~dodstp proc~djaccd->proc~derstep proc~djaccd->proc~dhstep proc~djacfd->proc~derstep proc~djacfd->proc~dhstep proc~djckc->proc~dpvb proc~djckc->proc~dpvd proc~djckf djckf proc~djckc->proc~djckf proc~djckz->proc~dpvb proc~djckz->proc~dpvd proc~dodpc1->proc~dhstep proc~dppt dppt proc~dodpc3->proc~dppt proc~dodstp->interface~dnrm2 proc~dodstp->proc~dwght proc~dodstp->proc~dfctr dchex dchex proc~dodstp->dchex dqrdc dqrdc proc~dodstp->dqrdc dqrsl dqrsl proc~dodstp->dqrsl dtrco dtrco proc~dodstp->dtrco dtrsl dtrsl proc~dodstp->dtrsl interface~drot drot proc~dodstp->interface~drot interface~drotg drotg proc~dodstp->interface~drotg interface~idamax idamax proc~dodstp->interface~idamax proc~desubi desubi proc~dodstp->proc~desubi proc~dsolve dsolve proc~dodstp->proc~dsolve proc~dvevtr dvevtr proc~dodstp->proc~dvevtr proc~djckf->proc~dpvb proc~djckf->proc~dpvd proc~dppnml dppnml proc~dppt->proc~dppnml proc~dsolve->interface~ddot interface~daxpy daxpy proc~dsolve->interface~daxpy proc~dvevtr->proc~dsolve

Variables

Type Attributes Name Initial
real(kind=wp) :: beta(3) = [1.1E-0_wp, 3.3E+0_wp, 8.7_wp]
integer :: lunrpt

Source Code

program example4
!! Default ODR job, with parameter bounds.
!!   This sample problem comes from Zwolak et al. 2001 (High Performance Computing
!! Symposium, "Estimating rate constants in cell cycle models"). The call to
!! ODRPACK95 is modified from the call the authors make to ODRPACK. This is
!! done to illustrate the need for bounds. The authors could just have easily
!! used the call statement here to solve their problem.
!!   Curious users are encouraged to remove the bounds in the call statement,
!! run the code, and compare the results to the current call statement.
   use odrpack_kinds, only: wp
   use odrpack, only: odr
   use example4_model, only: fcn
   implicit none

   real(kind=wp) :: beta(3) = [1.1E-0_wp, 3.3E+0_wp, 8.7_wp]
   integer :: lunrpt
   ! INTEGER :: I
   ! REAL (KIND=wp) :: C, M, TOUT

   open (newunit=lunrpt, file="./example/report4.dat")

   call odr(fcn, &
            n=5, m=1, np=3, nq=1, &
            beta=beta, &
            y=reshape([55.0_wp, 45.0_wp, 40.0_wp, 30.0_wp, 20.0_wp], [5, 1]), &
            x=reshape([0.15_wp, 0.20_wp, 0.25_wp, 0.30_wp, 0.50_wp], [5, 1]), &
            lower=[0.0_wp, 0.0_wp, 0.0_wp], &
            upper=[1000.0_wp, 1000.0_wp, 1000.0_wp], &
            iprint=2122, &
            lunrpt=lunrpt, &
            maxit=20)

   close (9)

   ! The following code will reproduce the plot in Figure 2 of Zwolak et al. 2001.
   !    DO I = 0, 100
   !       C = 0.05 + (0.7 - 0.05)*I/100
   !       TOUT = 1440.0D0
   !       !CALL MPF(M,C,1.1D-10,3.3D-3,8.7D0,0.0D0,TOUT,C/2)
   !       CALL MPF(M, C, 1.15395968E-02_wp, 2.61676386E-03_wp, &
   !                9.23138811E+00_wp, 0.0D0, TOUT, C/2)
   !       WRITE (*, *) C, TOUT
   !    END DO

end program example4