!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Structured Markov Chains Solver [ SMCSolver ] ! ! Dario Bini, Beatrice Meini, Sergio Steffe' ! ! bini@dm.unipi.it, meini@dm.unipi.it, steffe@dm.unipi.it ! ! Dipartimento di Matematica "Leonida Tonelli" ! ! Largo Pontecorvo 5 ! ! 56127 Pisa ! ! Italy ! ! Version 1.3 - March 2007 ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! ! ! module smc_int interface subroutine drft(am1, a0, a1, v, drift) integer, parameter :: dp=kind(0.d0) real(dp),intent(in),dimension(:,:) :: am1, a0, a1 real(dp), intent(out) :: drift real(dp),dimension(:),intent(out) :: v end subroutine drft end interface interface subroutine gth(a,p) integer, parameter :: dp=kind(0.d0) real(dp),dimension(:,:),intent(in) :: a real(dp),dimension(:),intent(out) :: p end subroutine gth end interface interface subroutine crqbd(Am1, A0, A1, G, R, U, doshift, dogth, drift, nerror, maxit) integer, parameter :: dp=kind(0.d0) real(dp),dimension(:,:),intent(in) :: Am1, A0, A1 real(dp),dimension(:,:),intent(out) :: G real(dp),dimension(:,:),intent(out),optional :: R, U real(dp),intent(out),optional ::drift logical, optional,intent(in) :: doshift, dogth integer,optional,intent(in) :: nerror integer, optional,intent(in) :: maxit end subroutine crqbd end interface interface subroutine lrqbd(Am1, A0, A1, G, R, U, doshift, dogth, drift, nerror, maxit) integer, parameter :: dp=kind(0.d0) real(dp),dimension(:,:),intent(in) :: Am1, A0, A1 real(dp),dimension(:,:),intent(out),optional :: G, R, U real(dp),intent(out),optional :: drift logical, optional,intent(in) :: doshift, dogth integer,optional,intent(in) :: nerror integer, optional,intent(in) :: maxit end subroutine lrqbd end interface interface subroutine bgth(v,a,b) integer, parameter :: dp=kind(0.d0) real(dp),intent(in), dimension(:,:) :: a real(dp),intent(inout), dimension(:,:) :: b real(dp),intent(inout),dimension(:) :: v end subroutine bgth end interface interface SUBROUTINE shift(a,a1,drift) integer, parameter :: dp=kind(0.d0) REAL(dp), DIMENSION(:,:,:):: a REAL(dp), DIMENSION(:,:) :: a1 real(dp) :: drift end SUBROUTINE shift end interface interface subroutine gm1tomg1(A, AA, v) integer,parameter::dp=kind(0.d0) real(dp),dimension(:,:,:):: A,AA real(dp),dimension(:):: v end subroutine gm1tomg1 end interface interface subroutine rresidual(A,r,norm1) integer,parameter::dp=kind(0.d0) real(dp),dimension(:,:,:):: A real(dp),dimension(:,:):: R real(dp)::norm1 end subroutine rresidual end interface interface subroutine gtou(A0,A1,G,U) integer,parameter:: dp=kind(0.d0) real(dp),dimension(:,:) :: a0,a1,g,u end subroutine gtou end interface interface subroutine utor(A1,u,r) integer,parameter:: dp=kind(0.d0) real(dp),dimension(:,:) :: a1,u,r end subroutine utor end interface interface subroutine qbdrres(Am1,A0,A1,r,res) integer,parameter:: dp=kind(0.d0) real(dp),dimension(:,:) :: am1,a0,a1,r real(dp)::res end subroutine qbdrres end interface interface subroutine qbdures(Am1,A0,A1,u,res) integer,parameter:: dp=kind(0.d0) real(dp),dimension(:,:) :: am1,a0,a1,u real(dp)::res end subroutine qbdures end interface end module smc_int module smc_tools ! Default RC simple shift to 0 if drift<0 ! to infty if drift>0 , ! option no shift ! main_subroutine: se compare solo G calcola solo G ! G e R calcola G e R ! U ! pi ! output ulteriore: drift, residuo ! raffinamento con Newton se residuo e' grande (opzionale?) ! (da valutare con esperimenti) integer,parameter :: dp=kind(0.d0) integer :: info=0 end module smc_tools