N.I.M.R.O.D.  

lineSearch.f90

Go to the documentation of this file.
00001 !------------------------------------------------------------------------------
00002 ! N.I.M.R.O.D. - Normal approximation Inference in Models with Random
00003 !                effects based on Ordinary Differential equations
00004 !------------------------------------------------------------------------------
00005 !
00006 ! VERSION : 1.0
00007 !
00008 ! MODULE: lineSearch.f90 function SEARPAS
00041 
00042 subroutine searpas(vw,step,b,bh,m,delta,fim,namefunc)
00043     !
00044     !  MINIMISATION UNIDIMENSIONNELLE
00045     !
00046     use mpimod
00047     use debugage
00048     implicit none
00049 
00050     integer,intent(in)::m
00051     double precision,dimension(m),intent(in)::b
00052     double precision,intent(inout)::vw
00053     double precision,dimension(m),intent(in)::bh,delta
00054     double precision,intent(inout)::fim,step
00055     double precision::vlw,vlw1,vlw2,vlw3,vm,fi1,fi2,fi3
00056     integer::i
00057     double precision,external::namefunc
00058 
00059     vlw1=dlog(vw)
00060     vlw2=vlw1+step
00061     call valfpa(vlw1,fi1,b,bh,m,delta,namefunc)
00062     call valfpa(vlw2,fi2,b,bh,m,delta,namefunc)
00063     if((numproc.EQ.0).and.(debugSEARPAS))print*,'vlw1',dexp(vlw1),fi1
00064     if((numproc.EQ.0).and.(debugSEARPAS))print*,'vlw2',dexp(vlw2),fi2
00065 
00066     if(fi2.ge.fi1) then
00067         vlw3=vlw2
00068         vlw2=vlw1
00069         fi3=fi2
00070         fi2=fi1
00071         step=-step
00072 
00073         vlw1=vlw2+step
00074 
00075         call valfpa(vlw1,fi1,b,bh,m,delta,namefunc)
00076         if((numproc.EQ.0).and.(debugSEARPAS))print*,'vlw3',dexp(vlw1),fi1
00077         ! COMMENTED 01/01/2013 if(fi1.gt.fi2) goto 50
00078     else
00079         vlw=vlw1
00080         vlw1=vlw2
00081         vlw2=vlw
00082         fim=fi1
00083         fi1=fi2
00084         fi2=fim
00085     end if
00086 
00087     do i=1,40
00088         vlw3=vlw2
00089         vlw2=vlw1
00090         fi3=fi2
00091         fi2=fi1
00092 
00093         vlw1=vlw2+step
00094         call valfpa(vlw1,fi1,b,bh,m,delta,namefunc)
00095         if((numproc.EQ.0).and.(debugSEARPAS))print*,'vlw4',dexp(vlw1),fi1
00096         if(fi1.gt.fi2) goto 50
00097         if(fi1.eq.fi2) then
00098             fim=fi2
00099             vm=vlw2
00100             goto 100
00101         end if
00102     end do
00103 !
00104 !  PHASE 2 APPROXIMATION PAR QUADRIQUE
00105 !
00106 50 continue
00107    !
00108    !  CALCUL MINIMUM QUADRIQUE
00109    !
00110    vm=vlw2-step*(fi1-fi3)/(2.d0*(fi1-2.d0*fi2+fi3))
00111    call valfpa(vm,fim,b,bh,m,delta,namefunc)
00112    if((numproc.EQ.0).and.(debugSEARPAS))print*,'vlwend',dexp(vm),fim
00113    if(fim.le.fi2) goto 100
00114    vm=vlw2
00115    fim=fi2
00116 100 continue
00117     vw=dexp(vm)
00118 
00119     return
00120 
00121 end subroutine searpas
00122 
00123 !------------------------------------------------------------------------------
00124 ! N.I.M.R.O.D. - Normal approximation Inference in Models with Random
00125 !                effects based on Ordinary Differential equations
00126 !------------------------------------------------------------------------------
00127 !
00128 ! VERSION : 1.0
00129 !
00130 ! MODULE: lineSearch.f90 function valfpa
00158 subroutine valfpa(vw,fi,b,bk,m,delta,namefunc)
00159     use WorkingSharedValues
00160     use mpimod
00161     implicit none
00162 
00163     integer,intent(in)::m
00164     double precision,dimension(m),intent(in)::b,delta
00165     double precision,dimension(m),intent(out)::bk
00166     double precision,intent(out)::fi
00167     double precision::vw,z
00168     integer::i0,i
00169     double precision,external::namefunc
00170 
00171     z=0.d0
00172     i0=1
00173     do i=1,m
00174         bk(i)=b(i)+dexp(vw)*delta(i)
00175     end do
00176     fi=-namefunc(bk,m,i0,z,i0,z)
00177     if((numproc.EQ.0).and.(debugSEARPAS)) &
00178         print*,"--------------->",dexp(vw),fi
00179     return
00180 
00181 
00182 end subroutine valfpa