!The Full Polymorphic Package
!Copyright (C) Etienne Forest
! Taylor polymorphism at execution is based on an idea
! and C++ prototype developed by J. Bengtsson circa 1990
module polymorphic_taylor
use complex_taylor
implicit none
public
logical(lp),private,parameter::t=.true.,f=.false.
integer,private::NO,ND,ND2,NP,NDPT,NV !,lastmaster ! 2002.12.13
INTEGER,PRIVATE::NMAX_pol=100
! real(dp),PRIVATE::EPS=c_1d_6
! integer ent,exi
integer,private,parameter::m1=mmmmmm1,m2=mmmmmm2,m3=mmmmmm3,ms=mmmmmm4
integer,private,parameter:: m11=m1+ms*m1,m12=m1+ms*m2,m13=m1+ms*m3, &
m21=m2+ms*m1,m22=m2+ms*m2,m23=m2+ms*m3, &
m31=m3+ms*m1,m32=m3+ms*m2,m33=m3+ms*m3
logical(lp),private::old
private resetpoly,resetpolyn ,resetpoly0,resetpolyn0 ,allocpoly,allocpolyn,resetpoly_R,resetpoly_Rn
PRIVATE K_OPT,A_OPT,e_OPT,ke_OPT
private allocenv,allocenvn,resetenv,resetenvn
private equal,Dequaldacon,equaldacon,iequaldacon ,realEQUAL,singleequal
private taylorEQUAL,EQUALtaylor,complexreal_8
private real_8univ,univreal_8
PRIVATE real_8MAP,MAPREAL_8
private unaryADD,add,daddsc,dscadd,addsc,scadd,iaddsc,iscadd
private unarySUB,subs,dscsub,dsubsc,scsub,subsc,iscsub,isubsc
private mul,dmulsc,dscmul,mulsc,scmul,imulsc,iscmul
private div,ddivsc,dscdiv,divsc,scdiv,idivsc,iscdiv
private POW,POWR,POWR8
private dexpt,dcost,dcosdt,dsindt,dsint,dlogt,dsqrtt
private greaterthan,greatereq
private lessthan,dlessthansc,dsclessthan,lessthansc,sclessthan,ilessthansc,isclessthan
private igreatersc,iscgreater,dgreatersc,dscgreater,greatersc,scgreater
private dgreatereqsc,dscgreatereq,greatereqsc,scgreatereq,igreatereqsc,iscgreatereq
private abst,full_abst
private lesseq,dlesseqsc,dsclesseq,lesseqsc,sclesseq,ilesseqsc,isclesseq
private eq,deqsc,dsceq,eqsc,sceq,ieqsc,isceq
private neq,dneqsc,dscneq,neqsc,scneq,ineqsc,iscneq
PRIVATE GETCHARnd2,GETintnd2,getchar,GETint,GETORDER,CUTORDER
private real_8ENV,ENVreal_8 ,RADTAYLORENV_8, ENV_8RADTAYLOR ,beamENV_8,normal_p
private absoftdatan2dr,absoftdcosdr,absoftdsindr,absoftdtandr,absoftdatandr
!complex stuff
private datant,datanDt,datan2t,dasint,dacost,dtant,dtandt,datanht
private dcosht,dsinht,dtanht,SINX_XT,SINHX_XT,polymorpht
! PRIVATE EQUAL1D,EQUAL2D
! end complex stuff
private printpoly,printdouble,printsingle,dmulmapconcat
private line
character(120) line
!integer npol
!parameter (npol=20)
!INTEGER iasspol
!type(real_8) DUMMYpol,DUMpol(ndum)
!integer iasspol0
!private assp
PRIVATE SINH_HR
PRIVATE SIN_HR
! PROBE_8 STUFF
PRIVATE RADTAYLORprobe_8,beamprobe_8
real(dp) :: sinhx_x_min=1e-4_dp
real(dp) :: sinhx_x_minp=1.0_dp ! 1.e-9 !c_0_0001
INTERFACE assignment (=)
MODULE PROCEDURE EQUAL ! 2002.10.9
! MODULE PROCEDURE EQUAL1D ! 2004.7.10
! MODULE PROCEDURE EQUAL2D ! 2004.7.10
MODULE PROCEDURE complexreal_8
!zgoubi
MODULE PROCEDURE realEQUAL !
MODULE PROCEDURE singleequal
MODULE PROCEDURE taylorEQUAL !taylor=real_8
MODULE PROCEDURE EQUALtaylor !real_8= taylor ! here 2002.10.9
MODULE PROCEDURE Dequaldacon !
MODULE PROCEDURE equaldacon !
! For resetting
MODULE PROCEDURE Iequaldacon
MODULE PROCEDURE Iequaldaconn
! For universal_taylor
MODULE PROCEDURE real_8univ
MODULE PROCEDURE univreal_8 ! new sagan
! for damap
MODULE PROCEDURE MAPreal_8
MODULE PROCEDURE real_8MAP
! FOR RADIATION
MODULE PROCEDURE real_8ENV
MODULE PROCEDURE ENVreal_8
MODULE PROCEDURE RADTAYLORENV_8
MODULE PROCEDURE ENV_8RADTAYLOR
! allowing normalization of polymorphic arrays directly
MODULE PROCEDURE beamENV_8
MODULE PROCEDURE normal_p
MODULE PROCEDURE RADTAYLORprobe_8
MODULE PROCEDURE beamprobe_8
end INTERFACE
!@
!@
!@ |
!@
!@ + |
!@
!@
!@ REAL_8 |
!@
!@
!@
!@ Real(dp) |
!@
!@
!@ Real(sp) |
!@
!@
!@ Integer |
!@
!@
!@ |
!@
!@ REAL_8 |
!@
!@
!@
!@ add |
!@
!@
!@
!@ daddsc |
!@
!@
!@ ADDSC |
!@
!@
!@ IADDSC |
!@
!@
!@ |
!@
!@
!@ Real(dp) |
!@
!@
!@
!@ dscadd |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
!@ |
!@
!@ Real(sp) |
!@
!@
!@ SCADD |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
!@ |
!@
!@ Integer |
!@
!@
!@ ISCADD |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
INTERFACE OPERATOR (+)
MODULE PROCEDURE unaryADD !
MODULE PROCEDURE add !
MODULE PROCEDURE daddsc !
MODULE PROCEDURE dscadd !
MODULE PROCEDURE addsc !
MODULE PROCEDURE scadd !
MODULE PROCEDURE iaddsc !
MODULE PROCEDURE iscadd !
END INTERFACE
!@
!@
!@ |
!@
!@ - |
!@
!@
!@ REAL_8 |
!@
!@
!@
!@ Real(dp) |
!@
!@
!@ Real(sp) |
!@
!@
!@ Integer |
!@
!@
!@ |
!@
!@ REAL_8 |
!@
!@
!@
!@ SUBS |
!@
!@
!@
!@ dSUBsc |
!@
!@
!@ SUBSC |
!@
!@
!@ ISUBSC |
!@
!@
!@ |
!@
!@
!@ Real(dp) |
!@
!@
!@
!@ dscSUB |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
!@ |
!@
!@ Real(sp) |
!@
!@
!@ SCSUB |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
!@ |
!@
!@ Integer |
!@
!@
!@ ISCSUB |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
INTERFACE OPERATOR (-)
MODULE PROCEDURE unarySUB !
MODULE PROCEDURE subs !
MODULE PROCEDURE dscsub !
MODULE PROCEDURE dsubsc !
MODULE PROCEDURE subsc !
MODULE PROCEDURE scsub !
MODULE PROCEDURE isubsc !
MODULE PROCEDURE iscsub !
END INTERFACE
!@
!@
!@ |
!@
!@ * |
!@
!@
!@ REAL_8 |
!@
!@
!@
!@ Real(dp) |
!@
!@
!@ Real(sp) |
!@
!@
!@ Integer |
!@
!@
!@ |
!@
!@ REAL_8 |
!@
!@
!@
!@ MUL |
!@
!@
!@
!@ dMULsc |
!@
!@
!@ MULSC |
!@
!@
!@ IMULSC |
!@
!@
!@ |
!@
!@
!@ Real(dp) |
!@
!@
!@
!@ dscMUL |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
!@ |
!@
!@ Real(sp) |
!@
!@
!@ SCMUL |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
!@ |
!@
!@ Integer |
!@
!@
!@ ISCMUL |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
INTERFACE OPERATOR (*)
MODULE PROCEDURE mul !
MODULE PROCEDURE dmulsc !
MODULE PROCEDURE dscmul !
MODULE PROCEDURE mulsc !
MODULE PROCEDURE scmul !
MODULE PROCEDURE imulsc !
MODULE PROCEDURE iscmul !
MODULE PROCEDURE dmulmapconcat ! concatenation real_8*damap
END INTERFACE
!@
!@
!@ |
!@
!@ / |
!@
!@
!@ REAL_8 |
!@
!@
!@
!@ Real(dp) |
!@
!@
!@ Real(sp) |
!@
!@
!@ Integer |
!@
!@
!@ |
!@
!@ REAL_8 |
!@
!@
!@
!@ div |
!@
!@
!@
!@ dDIVsc |
!@
!@ DIVSC |
!@
!@ IDIVSC |
!@
!@
!@ |
!@
!@
!@ Real(dp) |
!@
!@
!@
!@ dscDIV |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
!@ |
!@
!@ Real(sp) |
!@
!@ SCDIV |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
!@ |
!@
!@ Integer |
!@
!@ ISCDIV |
!@
!@ F90 |
!@
!@ F90 |
!@
!@ F90 |
!@
!@
INTERFACE OPERATOR (/)
MODULE PROCEDURE div !
MODULE PROCEDURE ddivsc !
MODULE PROCEDURE dscdiv !
MODULE PROCEDURE divsc !
MODULE PROCEDURE scdiv !
MODULE PROCEDURE idivsc !
MODULE PROCEDURE iscdiv !
END INTERFACE
INTERFACE OPERATOR (**)
MODULE PROCEDURE POW !
MODULE PROCEDURE POWR !
MODULE PROCEDURE POWR8 !
END INTERFACE
! Logical Operators
INTERFACE OPERATOR (<)
MODULE PROCEDURE lessthan
MODULE PROCEDURE dlessthansc
MODULE PROCEDURE dsclessthan
MODULE PROCEDURE lessthansc
MODULE PROCEDURE sclessthan
MODULE PROCEDURE ilessthansc
MODULE PROCEDURE isclessthan
END INTERFACE
INTERFACE OPERATOR (>)
MODULE PROCEDURE greaterthan
MODULE PROCEDURE dgreatersc
MODULE PROCEDURE dscgreater
MODULE PROCEDURE greatersc
MODULE PROCEDURE scgreater
MODULE PROCEDURE igreatersc
MODULE PROCEDURE iscgreater
END INTERFACE
INTERFACE OPERATOR (>=)
MODULE PROCEDURE greatereq
MODULE PROCEDURE dgreatereqsc
MODULE PROCEDURE dscgreatereq
MODULE PROCEDURE greatereqsc
MODULE PROCEDURE scgreatereq
MODULE PROCEDURE igreatereqsc
MODULE PROCEDURE iscgreatereq
END INTERFACE
INTERFACE OPERATOR (<=)
MODULE PROCEDURE lesseq
MODULE PROCEDURE dlesseqsc
MODULE PROCEDURE dsclesseq
MODULE PROCEDURE lesseqsc
MODULE PROCEDURE sclesseq
MODULE PROCEDURE ilesseqsc
MODULE PROCEDURE isclesseq
END INTERFACE
INTERFACE OPERATOR (==)
MODULE PROCEDURE eq
MODULE PROCEDURE deqsc
MODULE PROCEDURE dsceq
MODULE PROCEDURE eqsc
MODULE PROCEDURE sceq
MODULE PROCEDURE ieqsc
MODULE PROCEDURE isceq
END INTERFACE
INTERFACE OPERATOR (/=)
MODULE PROCEDURE neq
MODULE PROCEDURE dneqsc
MODULE PROCEDURE dscneq
MODULE PROCEDURE neqsc
MODULE PROCEDURE scneq
MODULE PROCEDURE ineqsc
MODULE PROCEDURE iscneq
END INTERFACE
! New Operators
INTERFACE OPERATOR (.SUB.)
MODULE PROCEDURE GETORDER
MODULE PROCEDURE getchar
MODULE PROCEDURE GETint
END INTERFACE
INTERFACE OPERATOR (.CUT.)
MODULE PROCEDURE CUTORDER
END INTERFACE
INTERFACE OPERATOR (.PAR.)
MODULE PROCEDURE getcharnd2
MODULE PROCEDURE GETintnd2
END INTERFACE
! Intrinsic function (some are Fortran extension not always supported)
INTERFACE exp
MODULE PROCEDURE dexpt
END INTERFACE
INTERFACE dexp
MODULE PROCEDURE dexpt
END INTERFACE
INTERFACE cexp
MODULE PROCEDURE dexpt
END INTERFACE
INTERFACE cdexp
MODULE PROCEDURE dexpt
END INTERFACE
INTERFACE cos
MODULE PROCEDURE dcost
END INTERFACE
INTERFACE ccos
MODULE PROCEDURE dcost
END INTERFACE
INTERFACE cdcos
MODULE PROCEDURE dcost
END INTERFACE
INTERFACE dcos
MODULE PROCEDURE dcost
END INTERFACE
INTERFACE atanh
MODULE PROCEDURE datanht
END INTERFACE
INTERFACE tan
MODULE PROCEDURE dtant
END INTERFACE
INTERFACE dtan
MODULE PROCEDURE dtant
END INTERFACE
INTERFACE tand
MODULE PROCEDURE dtandt
! MODULE PROCEDURE absoftdtandr ! for non PC
END INTERFACE
INTERFACE dtand
MODULE PROCEDURE dtandt
! MODULE PROCEDURE absoftdtandr ! for non PC absoft
END INTERFACE
INTERFACE cosd
MODULE PROCEDURE dcosdt
! MODULE PROCEDURE absoftdcosdr ! for non PC
END INTERFACE
INTERFACE dcosd
MODULE PROCEDURE dcosdt
! MODULE PROCEDURE absoftdcosdr ! for non PC absoft
END INTERFACE
INTERFACE sin
MODULE PROCEDURE dsint
END INTERFACE
INTERFACE cdsin
MODULE PROCEDURE dsint
END INTERFACE
INTERFACE ccsin
MODULE PROCEDURE dsint
END INTERFACE
INTERFACE dsin
MODULE PROCEDURE dsint
END INTERFACE
INTERFACE sind
MODULE PROCEDURE dsindt
! MODULE PROCEDURE absoftdsindr ! for non PC
END INTERFACE
INTERFACE dsind
MODULE PROCEDURE dsindt
! MODULE PROCEDURE absoftdsindr ! for non PC absoft
END INTERFACE
INTERFACE log
MODULE PROCEDURE dlogt
END INTERFACE
INTERFACE dlog
MODULE PROCEDURE dlogt
END INTERFACE
INTERFACE cdlog
MODULE PROCEDURE dlogt
END INTERFACE
INTERFACE clog
MODULE PROCEDURE dlogt
END INTERFACE
INTERFACE sqrt
MODULE PROCEDURE dsqrtt
END INTERFACE
INTERFACE dsqrt
MODULE PROCEDURE dsqrtt
END INTERFACE
INTERFACE abs
MODULE PROCEDURE abst
END INTERFACE
INTERFACE dabs
MODULE PROCEDURE abst
END INTERFACE
! complex stuff
INTERFACE atan
MODULE PROCEDURE datant
END INTERFACE
INTERFACE datan
MODULE PROCEDURE datant
END INTERFACE
INTERFACE atan2
MODULE PROCEDURE datan2t
END INTERFACE
INTERFACE datan2
MODULE PROCEDURE datan2t
END INTERFACE
INTERFACE atan2d
MODULE PROCEDURE datan2dt
! MODULE PROCEDURE absoftdatan2dr ! for non PC
END INTERFACE
INTERFACE datan2d
MODULE PROCEDURE datan2dt
! MODULE PROCEDURE absoftdatan2dr ! for non PC absoft
END INTERFACE
INTERFACE atand
MODULE PROCEDURE datandt
! MODULE PROCEDURE absoftdatandr ! for non PC
END INTERFACE
INTERFACE datand
MODULE PROCEDURE datandt
! MODULE PROCEDURE absoftdatandr ! for non PC absoft
END INTERFACE
INTERFACE asin
MODULE PROCEDURE dasint
END INTERFACE
INTERFACE dasin
MODULE PROCEDURE dasint
END INTERFACE
INTERFACE acos
MODULE PROCEDURE dacost
END INTERFACE
INTERFACE dacos
MODULE PROCEDURE dacost
END INTERFACE
INTERFACE cosh
MODULE PROCEDURE dcosht
END INTERFACE
INTERFACE dcosh
MODULE PROCEDURE dcosht
END INTERFACE
INTERFACE sinh
MODULE PROCEDURE dsinht
END INTERFACE
INTERFACE dsinh
MODULE PROCEDURE dsinht
END INTERFACE
INTERFACE tanh
MODULE PROCEDURE dtanht
END INTERFACE
INTERFACE dtanh
MODULE PROCEDURE dtanht
END INTERFACE
! end Intrinsic function
! Non Intrinsic function
INTERFACE full_abs
MODULE PROCEDURE full_abst
END INTERFACE
INTERFACE morph
MODULE PROCEDURE polymorpht
END INTERFACE
INTERFACE SINHX_X
MODULE PROCEDURE SINH_HR
MODULE PROCEDURE SINHX_XT
END INTERFACE
INTERFACE SINX_X
MODULE PROCEDURE SIN_HR
MODULE PROCEDURE SINX_XT
END INTERFACE
! i/o
INTERFACE daprint
MODULE PROCEDURE printpoly !
MODULE PROCEDURE printdouble
MODULE PROCEDURE printsingle
END INTERFACE
INTERFACE print
MODULE PROCEDURE printpoly !
MODULE PROCEDURE printdouble
MODULE PROCEDURE printsingle
END INTERFACE
! constructors and destructors
INTERFACE alloc
MODULE PROCEDURE A_OPT !allocpoly !
MODULE PROCEDURE allocpolyn !
!radiation
MODULE PROCEDURE e_opt !
! MODULE PROCEDURE allocenv !
MODULE PROCEDURE allocenvn !
END INTERFACE
INTERFACE kill
MODULE PROCEDURE K_OPT !resetpoly0 !
MODULE PROCEDURE Ke_OPT !resetpoly0 !
MODULE PROCEDURE resetpolyn0 !
MODULE PROCEDURE resetpoly_R
MODULE PROCEDURE resetpoly_RN
!radiation
! MODULE PROCEDURE resetenv !
MODULE PROCEDURE resetenvn !
END INTERFACE
INTERFACE reset
MODULE PROCEDURE resetpoly !
MODULE PROCEDURE resetpolyn !
!radiation
MODULE PROCEDURE resetenv !
MODULE PROCEDURE resetenvn !
END INTERFACE
! Managing
INTERFACE ass
MODULE PROCEDURE assp
END INTERFACE
contains
subroutine make_it_knob(k,i,s)
implicit none
TYPE (real_8), intent(inout) :: k
real(dp), optional :: s
integer, intent(in) :: i
k%s=1.0_dp
if(present(s)) k%s=s
k%i=i
k%kind=3
end subroutine make_it_knob
subroutine kill_knob(k)
implicit none
TYPE (real_8), intent(inout) :: k
k%s=1.0_dp
k%i=0
k%kind=1
end subroutine kill_knob
FUNCTION polymorpht( S1 )
implicit none
TYPE (real_8) polymorpht
TYPE (taylor), INTENT (IN) :: S1
integer localmaster
localmaster=master
call ass(polymorpht)
polymorpht%t= s1
master=localmaster
END FUNCTION polymorpht
FUNCTION GETCHARnd2( S1, S2 )
implicit none
TYPE (real_8) GETCHARnd2
TYPE (real_8), INTENT (IN) :: S1
CHARACTER(*) , INTENT (IN) :: S2
integer localmaster
type(taylor) t
localmaster=master
call ass(GETCHARnd2)
call alloc(t)
t=s1
t=t.par.s2
GETCHARnd2%t=t
! GETCHARnd2%kind=m2
! if(s1%kind==m2) then
! GETCHARnd2%t=s1%t.par.s2
! else
! GETCHARnd2%t=zero
! endif
call kill(t)
master=localmaster
END FUNCTION GETCHARnd2
FUNCTION GETintnd2( S1, S2 )
implicit none
TYPE (real_8) GETintnd2
TYPE (real_8), INTENT (IN) :: S1
integer, INTENT (IN) :: S2(:)
integer localmaster
type(taylor) t
localmaster=master
call ass(GETintnd2)
call alloc(t)
t=s1
! if(s1%kind==m2) then
! GETintnd2%t=s1%t.par.s2
! else
! GETintnd2%kind=m2
! GETintnd2%t=zero
!
! endif
t=t.par.s2
GETintnd2%t=t
call kill(t)
master=localmaster
END FUNCTION GETintnd2
FUNCTION GETchar( S1, S2 )
implicit none
real(dp) GETchar
TYPE (real_8), INTENT (IN) :: S1
CHARACTER(*) , INTENT (IN) :: S2
! integer localmaster
GETchar=0.0_dp
if(s1%kind==m2) then
! GETchar%t=s1%t.sub.s2 ! OLD
GETchar=s1%t.sub.s2 ! CHANGE
endif
END FUNCTION GETchar
FUNCTION GETint( S1, S2 )
implicit none
real(dp) GETint
TYPE (real_8), INTENT (IN) :: S1
integer , INTENT (IN) :: S2(:)
integer i
GETint=0.0_dp
if(s1%kind==m2) then
GETint=s1%t.sub.s2 ! CHANGE
elseif(s1%kind==m1) then
!zgoubi
GETint=s1
do i=1,size(s2)
if(S2(i)/=0) then
GETint=0.0_dp
exit
endif
enddo
endif
END FUNCTION GETint
FUNCTION GETORDER( S1, S2 )
implicit none
TYPE (real_8) GETORDER
TYPE (real_8), INTENT (IN) :: S1
integer , INTENT (IN) :: S2
integer localmaster
localmaster=master
call ass(GETORDER)
if(s1%kind==m2) then
GETORDER%t=s1%t.sub.s2
else
GETORDER%kind=m1
GETORDER%r=0.0_dp
if(s2==0) GETORDER%r=s1%r
endif
master=localmaster
END FUNCTION GETORDER
FUNCTION CUTORDER( S1, S2 )
implicit none
TYPE (real_8) CUTORDER
TYPE (real_8), INTENT (IN) :: S1
integer , INTENT (IN) :: S2
integer localmaster
localmaster=master
call ass(CUTORDER)
CUTORDER=0.0_dp
if(s1%kind==m2) then
cutorder%kind=m2
CUTORDER%t=s1%t.CUT.s2
elseif(s1%kind==m1) then
if(s2>=1) CUTORDER%r=s1%r
cutorder%kind=m1
endif
master=localmaster
END FUNCTION CUTORDER
FUNCTION greatereq( S1, S2 )
implicit none
logical(lp) greatereq
TYPE (real_8), INTENT (IN) :: S1,S2
greatereq=.FALSE.
select case(s1%kind+ms*s2%kind)
case(m11)
IF(S1%r>=S2%r) greatereq=.TRUE.
case(m12,m21,m22)
select case(s1%kind+ms*s2%kind)
case(m21)
IF((S1%t.SUB.'0')>=S2%r) greatereq=.TRUE.
case(m12)
IF(S1%r>=(S2%t.SUB.'0')) greatereq=.TRUE.
case(m22)
IF((S1%t.SUB.'0')>=(S2%t.SUB.'0')) greatereq=.TRUE.
end select
case(m13,m31,m32,m23,m33)
select case(s1%kind+ms*s2%kind)
case(m31,m13,m33)
IF(S1%r>=S2%r) greatereq=.TRUE.
case(m32)
IF(S1%r>=(S2%t.SUB.'0')) greatereq=.TRUE.
case(m23)
IF((S1%t.SUB.'0')>=S2%r) greatereq=.TRUE.
end select
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in greatereq "
w_p%c(2)= "s1%kind ,s2%kind "
w_p=(/s1%kind ,s2%kind/)
! call !write_e(0)
end select
END FUNCTION greatereq
FUNCTION dgreatereqsc( S1, S2 )
implicit none
logical(lp) dgreatereqsc
TYPE (real_8), INTENT (IN) :: S1
real(dp), INTENT (IN) :: S2
dgreatereqsc=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S1%r>=S2) dgreatereqsc=.TRUE.
case(m2)
IF((S1%t.SUB.'0')>=S2) dgreatereqsc=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(1((1X,i4)))'
w_p%c(1)= " trouble in dgreatereqsc "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind/)
! call !write_e(0)
end select
END FUNCTION dgreatereqsc
FUNCTION dscgreatereq( S2,S1 )
implicit none
logical(lp) dscgreatereq
TYPE (real_8), INTENT (IN) :: S1
real(dp), INTENT (IN) :: S2
dscgreatereq=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S2>=S1%r) dscgreatereq=.TRUE.
case(m2)
IF(S2>=(S1%t.SUB.'0')) dscgreatereq=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(1((1X,i4)))'
w_p%c(1)= " trouble in dscgreatereq "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind/)
! call !write_e(0)
end select
END FUNCTION dscgreatereq
FUNCTION greatereqsc( S1, S2 )
implicit none
logical(lp) greatereqsc
TYPE (real_8), INTENT (IN) :: S1
real(sp), INTENT (IN) :: S2
if(real_warning) call real_stop
greatereqsc=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S1%r>=S2) greatereqsc=.TRUE.
case(m2)
IF((S1%t.SUB.'0')>=S2) greatereqsc=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(1((1X,i4)))'
w_p%c(1)= " trouble in greatereqsc "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind/)
! call !write_e(0)
end select
END FUNCTION greatereqsc
FUNCTION scgreatereq( S2,S1 )
implicit none
logical(lp) scgreatereq
TYPE (real_8), INTENT (IN) :: S1
real(sp), INTENT (IN) :: S2
if(real_warning) call real_stop
scgreatereq=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S2>=S1%r) scgreatereq=.TRUE.
case(m2)
IF(S2>=(S1%t.SUB.'0')) scgreatereq=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(1((1X,i4)))'
w_p%c(1)= " trouble in scgreatereq "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind/)
! call !write_e(0)
end select
END FUNCTION scgreatereq
FUNCTION igreatereqsc( S1, S2 )
implicit none
logical(lp) igreatereqsc
TYPE (real_8), INTENT (IN) :: S1
integer, INTENT (IN) :: S2
igreatereqsc=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S1%r>=S2) igreatereqsc=.TRUE.
case(m2)
IF((S1%t.SUB.'0')>=S2) igreatereqsc=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(1((1X,i4)))'
w_p%c(1)= " trouble in igreatereqsc "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind/)
! call !write_e(0)
end select
END FUNCTION igreatereqsc
FUNCTION iscgreatereq( S2,S1 )
implicit none
logical(lp) iscgreatereq
TYPE (real_8), INTENT (IN) :: S1
integer, INTENT (IN) :: S2
iscgreatereq=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S2>=S1%r) iscgreatereq=.TRUE.
case(m2)
IF(S2>=(S1%t.SUB.'0')) iscgreatereq=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(1((1X,i4)))'
w_p%c(1)= " trouble in iscgreatereq "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind/)
! call !write_e(0)
end select
END FUNCTION iscgreatereq
FUNCTION greaterthan( S1, S2 )
implicit none
logical(lp) greaterthan
TYPE (real_8), INTENT (IN) :: S1,S2
greaterthan=.FALSE.
select case(s1%kind+ms*s2%kind)
case(m11)
IF(S1%r>S2%r) greaterthan=.TRUE.
case(m12,m21,m22)
select case(s1%kind+ms*s2%kind)
case(m21)
IF((S1%t.SUB.'0')>S2%r) greaterthan=.TRUE.
case(m12)
IF(S1%r>(S2%t.SUB.'0')) greaterthan=.TRUE.
case(m22)
IF((S1%t.SUB.'0')>(S2%t.SUB.'0')) greaterthan=.TRUE.
end select
case(m13,m31,m32,m23,m33)
select case(s1%kind+ms*s2%kind)
case(m31,m13,m33)
IF(S1%r>S2%r) greaterthan=.TRUE.
case(m32)
IF(S1%r>(S2%t.SUB.'0')) greaterthan=.TRUE.
case(m23)
IF((S1%t.SUB.'0')>S2%r) greaterthan=.TRUE.
end select
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in greaterthan "
w_p%c(2)= "s1%kind ,s2%kind "
w_p=(/s1%kind ,s2%kind/)
! call !write_e(0)
end select
END FUNCTION greaterthan
FUNCTION igreatersc( S1, S2 )
implicit none
logical(lp) igreatersc
TYPE (real_8), INTENT (IN) :: S1
integer, INTENT (IN) :: S2
igreatersc=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S1%r>S2) igreatersc=.TRUE.
case(m2)
IF((S1%t.SUB.'0')>S2) igreatersc=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(1((1X,i4)))'
w_p%c(1)= " trouble in igreatersc "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind/)
! call !write_e(0)
end select
END FUNCTION igreatersc
FUNCTION iscgreater( S2,S1 )
implicit none
logical(lp) iscgreater
TYPE (real_8), INTENT (IN) :: S1
integer, INTENT (IN) :: S2
iscgreater=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S2>S1%r) iscgreater=.TRUE.
case(m2)
IF(S2>(S1%t.SUB.'0')) iscgreater=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(1((1X,i4)))'
w_p%c(1)= " trouble in iscgreater "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind/)
! call !write_e(0)
end select
END FUNCTION iscgreater
FUNCTION dgreatersc( S1, S2 )
implicit none
logical(lp) dgreatersc
TYPE (real_8), INTENT (IN) :: S1
real(dp), INTENT (IN) :: S2
dgreatersc=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S1%r>S2) dgreatersc=.TRUE.
case(m2)
IF((S1%t.SUB.'0')>S2) dgreatersc=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(1((1X,i4)))'
w_p%c(1)= " trouble in dgreatersc "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind/)
! call !write_e(0)
end select
END FUNCTION dgreatersc
FUNCTION dscgreater( S2,S1 )
implicit none
integer ipause, mypause
logical(lp) dscgreater
TYPE (real_8), INTENT (IN) :: S1
real(dp), INTENT (IN) :: S2
dscgreater=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S2>S1%r) dscgreater=.TRUE.
case(m2)
IF(S2>(S1%t.SUB.'0')) dscgreater=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(1((1X,i4)))'
w_p%c(1)= " trouble in dscgreater "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind/)
! call !write_e(0)
ipause=mypause(0)
end select
END FUNCTION dscgreater
FUNCTION greatersc( S1, S2 )
implicit none
logical(lp) greatersc
TYPE (real_8), INTENT (IN) :: S1
real(sp), INTENT (IN) :: S2
if(real_warning) call real_stop
greatersc=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S1%r>S2) greatersc=.TRUE.
case(m2)
IF((S1%t.SUB.'0')>S2) greatersc=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in greatersc "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END FUNCTION greatersc
FUNCTION scgreater( S2,S1 )
implicit none
logical(lp) scgreater
TYPE (real_8), INTENT (IN) :: S1
real(sp), INTENT (IN) :: S2
if(real_warning) call real_stop
scgreater=.FALSE.
select case(s1%kind)
case(m1,m3)
IF(S2>S1%r) scgreater=.TRUE.
case(m2)
IF(S2>(S1%t.SUB.'0')) scgreater=.TRUE.
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in scgreater "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END FUNCTION scgreater
FUNCTION lessthan( S1, S2 )
implicit none
logical(lp) lessthan
TYPE (real_8), INTENT (IN) :: S1,S2
lessthan=.FALSE.
select case(s1%kind+ms*s2%kind)
case(m11)
IF(S1%r0) then
write(i,*) s2%r," +",s2%s," *x_",int(s2%i)
else
write(i,*) s2%r
endif
if(s2%alloc) then
write(line,'(a41)') " weird Taylor part should be deallocated "
ipause=mypauses(0,line)
endif
end select
else
line= "Warning not defined in Printpoly (real polymorph)"
ipause=mypauses(1,line)
endif
END SUBROUTINE printpoly
SUBROUTINE printdouble(S2,i)
implicit none
real(dp),INTENT(INOUT)::S2
integer i
write(i,*) s2
END SUBROUTINE printdouble
SUBROUTINE printsingle(S2,i)
implicit none
real(sp),INTENT(INOUT)::S2
integer i
write(i,*) s2
END SUBROUTINE printsingle
SUBROUTINE resetpoly(S2)
implicit none
type (real_8),INTENT(INOUT)::S2
if(s2%alloc) call killtpsa(s2%t)
s2%alloc=f
s2%kind=0
s2%r=0.0_dp
!s2%s=one
!s2%i=0
END SUBROUTINE resetpoly
SUBROUTINE resetpolyn(S2,K)
implicit none
type (real_8),INTENT(INOUT),dimension(:)::S2
INTEGER,optional,INTENT(IN)::k
INTEGER J,i,N
if(present(k)) then
I=LBOUND(S2,DIM=1)
N=LBOUND(S2,DIM=1)+K-1
else
I=LBOUND(S2,DIM=1)
N=UBOUND(S2,DIM=1)
endif
DO J=I,N
call resetpoly(s2(j))
enddo
END SUBROUTINE resetpolyn
SUBROUTINE resetpoly_R(S2,FL) ! STAYS REAL
implicit none
type (real_8),INTENT(INOUT)::S2
logical(lp),INTENT(IN)::FL
if(s2%alloc) call killtpsa(s2%t)
s2%alloc=f
s2%kind=1
s2%r=0.0_dp
IF(.NOT.FL) THEN
s2%i=0
s2%s=1.0_dp
ENDIF
END SUBROUTINE resetpoly_R
SUBROUTINE resetpoly_RN(S2,FL,K)
implicit none
type (real_8),INTENT(INOUT),dimension(:)::S2
logical(lp),INTENT(IN)::FL
INTEGER,optional,INTENT(IN)::k
INTEGER J,i,N
if(present(k)) then
I=LBOUND(S2,DIM=1)
N=LBOUND(S2,DIM=1)+K-1
else
I=LBOUND(S2,DIM=1)
N=UBOUND(S2,DIM=1)
endif
DO J=I,N
call resetpoly_R(s2(j),FL)
enddo
END SUBROUTINE resetpoly_RN
SUBROUTINE resetpoly_R31(S2) ! STAYS REAL FOR PTC
implicit none
integer ipause, mypauses
type (real_8),INTENT(INOUT)::S2
if(s2%kind==3) then
if(s2%alloc) then
line= "Allocated in resetpoly_R31"
ipause=mypauses(2,line)
endif
s2%kind=1
s2%i=0
s2%s=1.0_dp
ENDIF
END SUBROUTINE resetpoly_R31
SUBROUTINE resetpoly_R31N(S2,K)
implicit none
type (real_8),INTENT(INOUT),dimension(:)::S2
INTEGER,optional,INTENT(IN)::k
INTEGER J,i,N
if(present(k)) then
I=LBOUND(S2,DIM=1)
N=LBOUND(S2,DIM=1)+K-1
else
I=LBOUND(S2,DIM=1)
N=UBOUND(S2,DIM=1)
endif
DO J=I,N
call resetpoly_R31(s2(j))
enddo
END SUBROUTINE resetpoly_R31N
SUBROUTINE resetpoly0(S2)
implicit none
type (real_8),INTENT(INOUT)::S2
if(s2%alloc) call killtpsa(s2%t)
s2%alloc=f
s2%kind=0
s2%r=0.0_dp
s2%i=0
s2%s=1.0_dp
END SUBROUTINE resetpoly0
SUBROUTINE K_OPT(S1,S2,s3,s4,s5,s6,s7,s8,s9,s10)
implicit none
type (real_8),INTENT(INout)::S1
type (real_8),optional, INTENT(INout):: S2,s3,s4,s5,s6,s7,s8,s9,s10
call resetpoly0(s1)
if(present(s2)) call resetpoly0(s2)
if(present(s3)) call resetpoly0(s3)
if(present(s4)) call resetpoly0(s4)
if(present(s5)) call resetpoly0(s5)
if(present(s6)) call resetpoly0(s6)
if(present(s7)) call resetpoly0(s7)
if(present(s8)) call resetpoly0(s8)
if(present(s9)) call resetpoly0(s9)
if(present(s10))call resetpoly0(s10)
END SUBROUTINE K_OPT
SUBROUTINE Ke_OPT(S1,S2,s3,s4,s5,s6,s7,s8,s9,s10)
implicit none
type (env_8),INTENT(INout)::S1
type (env_8),optional, INTENT(INout):: S2,s3,s4,s5,s6,s7,s8,s9,s10
call resetenv(s1)
if(present(s2)) call resetenv(s2)
if(present(s3)) call resetenv(s3)
if(present(s4)) call resetenv(s4)
if(present(s5)) call resetenv(s5)
if(present(s6)) call resetenv(s6)
if(present(s7)) call resetenv(s7)
if(present(s8)) call resetenv(s8)
if(present(s9)) call resetenv(s9)
if(present(s10))call resetenv(s10)
END SUBROUTINE Ke_OPT
! SUBROUTINE kill_c(k)
! implicit none
! integer ipause, mypauses
! integer k
! call dallsta(exi) ! Etienne checking routine
! if(exi/=ent) then
! write(line,'(3(1x,i8))') ent,exi,k
! ipause=mypauses(1999,line)
! k=k*10000
! endif
! END SUBROUTINE kill_c
SUBROUTINE resetpolyn0(S2,K)
implicit none
type (real_8),INTENT(INOUT),dimension(:)::S2
INTEGER,optional,INTENT(IN)::k
INTEGER J,i,N
if(present(k)) then
I=LBOUND(S2,DIM=1)
N=LBOUND(S2,DIM=1)+K-1
else
I=LBOUND(S2,DIM=1)
N=UBOUND(S2,DIM=1)
endif
DO J=I,N
call resetpoly0(s2(j))
enddo
END SUBROUTINE resetpolyn0
SUBROUTINE allocpoly(S2)
implicit none
type (real_8),INTENT(INOUT)::S2
!if(s2%alloc) call killtpsa(s2%t) ADDED ETIENNE
s2%alloc=f
s2%kind=1
s2%r=0.0_dp
s2%i=0
s2%g=0
s2%nb=0
s2%s=1.0_dp
END SUBROUTINE allocpoly
SUBROUTINE A_OPT(S1,S2,s3,s4,s5,s6,s7,s8,s9,s10)
implicit none
type (real_8),INTENT(INout)::S1
type (real_8),optional, INTENT(INout):: S2,s3,s4,s5,s6,s7,s8,s9,s10
call allocpoly(s1)
if(present(s2)) call allocpoly(s2)
if(present(s3)) call allocpoly(s3)
if(present(s4)) call allocpoly(s4)
if(present(s5)) call allocpoly(s5)
if(present(s6)) call allocpoly(s6)
if(present(s7)) call allocpoly(s7)
if(present(s8)) call allocpoly(s8)
if(present(s9)) call allocpoly(s9)
if(present(s10))call allocpoly(s10)
END SUBROUTINE A_opt
! SUBROUTINE alloc_c
! implicit none
! Etienne checking routine
! call dallsta(ent)
! END SUBROUTINE alloc_c
SUBROUTINE allocpolyn(S2,K)
implicit none
type (real_8),INTENT(INOUT),dimension(:)::S2
INTEGER,optional,INTENT(IN)::k
INTEGER J,i,N
if(present(k)) then
I=LBOUND(S2,DIM=1)
N=LBOUND(S2,DIM=1)+K-1
else
I=LBOUND(S2,DIM=1)
N=UBOUND(S2,DIM=1)
endif
DO J=I,N
call allocpoly(s2(j))
enddo
END SUBROUTINE allocpolyn
subroutine init_map_p(NO1,ND1,NP1,NDPT1,log)
implicit none
integer NO1,ND1,NP1,NDPT1
logical(lp) log
call init_map_c(NO1,ND1,NP1,NDPT1,log)
call set_in_polyp(log)
end subroutine init_map_p
subroutine init_tpsa_p(NO1,NP1,log)
implicit none
integer NO1,NP1
logical(lp) log
call init_tpsa_c(NO1,NP1,log)
call set_in_polyp(log)
end subroutine init_tpsa_p
subroutine set_in_polyp(log)
implicit none
logical(lp) log
integer iia(4),icoast(4)
call liepeek(iia,icoast)
old=log
NO=iia(1)
ND=iia(3)
ND2=iia(3)*2
NP=iia(2)-nd2
NDPT=icoast(4)
NV=iia(2)
! i_ =cmplx(zero,one,kind=dp)
end subroutine set_in_polyp
SUBROUTINE EQUAL2D(S2,S1)
implicit none
type (real_8),INTENT(inOUT)::S2(:,:)
type (real_8),INTENT(IN)::S1(:,:)
integer i,J,I1(2),I2(2)
I1(1)=LBOUND(S1,DIM=1)
I1(2)=LBOUND(S1,DIM=2)
I2(1)=LBOUND(S2,DIM=1)
I2(2)=LBOUND(S2,DIM=2)
do i=I1(1),UBOUND(S1,DIM=1)
do J=I1(2),UBOUND(S1,DIM=2)
S2(I-I1(1)+I2(1),J-I1(2)+I2(2))=S1(I,J)
ENDDO
ENDDO
end SUBROUTINE EQUAL2D
SUBROUTINE EQUAL1D(S2,S1)
implicit none
type (real_8),INTENT(inOUT)::S2(:)
type (real_8),INTENT(IN)::S1(:)
integer i,I1,I2
I1=LBOUND(S1,DIM=1)
I2=LBOUND(S2,DIM=1)
do i=I1,UBOUND(S1,DIM=1)
S2(I-I1+I2)=S1(I)
ENDDO
end SUBROUTINE EQUAL1D
SUBROUTINE EQUAL(S2,S1)
implicit none
integer ipause, mypauses
type (real_8),INTENT(inOUT)::S2
type (real_8),INTENT(IN)::S1
! integer localmaster
if(s1%kind==0) then
line= " You are putting kind=0 into something"
ipause=mypauses(10,line)
endif
if(s2%kind==3.and.(.not.setknob)) then
line= " You are putting something into a knob kind=3"
ipause=mypauses(11,line)
endif
if (s2%kind>0) then ! S2 exist
if(S2%kind==S1%kind) then
select case(S1%kind)
case(m1)
! localmaster=master ! added because of knob problems 2001.9.19
! call check_snake ! added because of knob problems 2001.9.19
! master=0
S2%R=S1%R
! master=localmaster ! added because of knob problems 2001.9.19
case(m2)
! localmaster=master
call check_snake
! 2002.12.26 master=0
S2%t=S1%t
! master=localmaster
case(m3)
s2%r=S1%r ! Knob stays a knob and real stays real 2002.10.9
! line= " You are putting kind=3 (TPSA) into another kind=3"
! ipause=mypauses(13,line)
end select
elseif(S2%kind>S1%kind ) then
if(S1%kind/=2) then
s2%r=S1%r
if(s2%kind/=3) s2%kind=1 !Knob stays a knob and real stays real 2002.10.9
else
s2%r=S1%t.sub.'0' ! setting a kind=3
endif
else
select case(s1%kind)
case(m2)
if(.not.s2%alloc) then
call alloc(s2%t)
s2%alloc=t
endif
s2%kind=2
! localmaster=master
call check_snake
! 2002.12.26 master=0
S2%t=S1%t
! master=localmaster
case(m3)
if(.not.s2%alloc) then
call alloc(s2%t)
s2%alloc=t
endif
s2%kind=2
! localmaster=master
call check_snake
! 2002.12.26 master=0
if(knob) then
call varfk1(S1)
S2%t=varf1
else
S2%R=S1%R
s2%kind=1
endif
! master=localmaster
end select
endif
else ! S2 does not exist
if(S1%kind==1) then ! what is s1
if(s2%i==0) then
S2%R=S1%R
s2%kind=1
elseif(s2%i>0.and.s2%i<=nv) then
call alloc(s2%t)
s2%t=(/S1%R,S2%S/).var.s2%i
! call var(s2%t,S1%R,S2%S,s2%i)
! s2%i=0
s2%kind=2
s2%alloc=t
else
write(6,*) "EQUAL IN m_POLYMORPH ",s2%i,S2%KIND,S2%ALLOC
WRITE(6,*) " I "
READ(5,*) s2%i
WRITE(6,*) SQRT(FLOAT(s2%i))
stop 777
endif
else
if(insane_PTC) then
if(.not.s2%alloc) then
call alloc(s2%t)
endif
! localmaster=master
call check_snake
! 2002.12.26 master=0
S2%t=S1%t
! master=localmaster
s2%kind=2
s2%alloc=t
else
w_p=0
w_p%nc=5
w_p%fc='(4(1X,A72,/),(1X,A72))'
write(w_p%c(1),'(A23,I4,A19)') " You are putting kind= ", s1%kind," (TPSA) in a kind=0"
w_p%c(2)= " We do not allow that anymore for safety reasons"
w_p%c(3)= " If you insist on it, modify real_polymorph and complex_polymorph"
w_p%c(4)= " at your own insane risk "
w_p%c(5)= " Etienne Forest/Frank Schmidt"
! call !write_e(778)
w_p%nc=sqrt(-float(w_p%nc))
endif
endif ! end of what is s1
endif ! S2 does not exist
END SUBROUTINE EQUAL
SUBROUTINE complexreal_8(S2,S1)
implicit none
complex(dp) ,INTENT(inOUT)::S2
type (real_8),INTENT(IN)::S1
! integer localmaster
select case(S1%kind)
case(m1)
S2=S1%R
case(m2)
! localmaster=master
call check_snake
! 2002.12.26 master=0
S2=S1%t.sub.'0'
! master=localmaster
case(m3)
! localmaster=master
call check_snake
! 2002.12.26 master=0
S2=S1%r
! master=localmaster
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in complexreal_8 "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END SUBROUTINE complexreal_8
SUBROUTINE realEQUAL(S2,S1)
implicit none
real(dp),INTENT(inOUT)::S2
type (real_8),INTENT(IN)::S1
! integer localmaster
select case(S1%kind)
case(m1)
S2=S1%R
case(m2)
! localmaster=master
call check_snake
! 2002.12.26 master=0
S2=S1%t.sub.'0'
! master=localmaster
case(m3)
! localmaster=master
call check_snake
! 2002.12.26 master=0
S2=S1%r
! master=localmaster
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in realEQUAL "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END SUBROUTINE realEQUAL
SUBROUTINE singleEQUAL(S2,S1)
implicit none
real(sp),INTENT(inOUT)::S2
type (real_8),INTENT(IN)::S1
! integer localmaster
select case(S1%kind)
case(m1)
S2=S1%R
case(m2)
! localmaster=master
call check_snake
! 2002.12.26 master=0
S2=S1%t.sub.'0'
! master=localmaster
case(m3)
! localmaster=master
call check_snake
! 2002.12.26 master=0
S2=S1%r
! master=localmaster
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in realEQUAL "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END SUBROUTINE singleEQUAL
SUBROUTINE taylorEQUAL(S2,S1)
implicit none
type (taylor),INTENT(inOUT)::S2
type (real_8),INTENT(IN)::S1
! integer localmaster
select case(S1%kind)
case(m1)
S2=S1%R
case(m2)
! localmaster=master
call check_snake
! 2002.12.26 master=0
S2=S1%t
! master=localmaster
case(m3)
! localmaster=master
call check_snake
! 2002.12.26 master=0
if(knob) then
call varfk1(S1)
S2=varf1
else
S2=S1%R
endif
! master=localmaster
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in taylorEQUAL "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END SUBROUTINE taylorEQUAL
SUBROUTINE EQUALtaylor(S2,S1)
implicit none
integer ipause, mypauses
type (real_8),INTENT(inOUT)::S2
type (taylor),INTENT(IN)::S1
! integer localmaster
IF(S2%KIND==3.and.(.not.setknob)) THEN ! 2002.10.9
line= "Forbidden in EQUALtaylor: s2 is a knob "
ipause=mypauses(0,line)
ENDIF
! localmaster=master
call check_snake
! 2002.12.26 master=0
if(s2%kind/=3) then
if(.not.s2%alloc) then
call alloc(s2%t)
s2%alloc=t
endif
S2%t=S1
s2%KIND=2
else
s2%r=S1.sub.'0' ! 2002.10.9
endif
! master=localmaster
END SUBROUTINE EQUALtaylor
SUBROUTINE real_8univ(S2,S1)
implicit none
integer ipause, mypauses
type (real_8),INTENT(inOUT)::S2
type (universal_taylor),INTENT(IN)::S1
! integer localmaster
IF(S2%KIND==M3) THEN
line= "Forbidden in real_8univ: s2 is a knob "
ipause=mypauses(0,line)
ENDIF
! localmaster=master
! call check_snake
! 2002.12.26 master=0
if(.not.s2%alloc) then
call alloc(s2%t)
s2%alloc=t
endif
S2%t=S1
s2%KIND=2
! master=localmaster
END SUBROUTINE real_8univ
SUBROUTINE univreal_8(S2,S1) ! new sagan
implicit none
type (universal_taylor),INTENT(inOUT)::S2
type (real_8),INTENT(IN)::S1
! integer localmaster
select case(S1%kind)
case(m1)
S2=S1%R
case(m2)
! localmaster=master
call check_snake
! 2002.12.26 master=0
S2=S1%t
! master=localmaster
case(m3)
! localmaster=master
call check_snake
! 2002.12.26 master=0
if(knob) then
call varfk1(S1)
S2=varf1
else
S2=S1%R
endif
! master=localmaster
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in univreal_8 "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END SUBROUTINE univreal_8
SUBROUTINE Dequaldacon(S2,R1)
implicit none
integer ipause, mypauses
type (real_8),INTENT(inOUT)::S2
real(dp),INTENT(IN)::R1
IF(S2%KIND==M3) THEN
if(setknob) then
s2%r=r1
return
else
line= "Forbidden in Dequaldacon: s2 is a knob "
ipause=mypauses(0,line)
endif
ENDIF
!localmaster=master
! call check_snake
! master=0
if (s2%kind/=0) then ! S2 exist
if(S2%kind==1.or.s2%kind==3) then
S2%R=R1
else
! S2%t=R1
s2%r=r1
s2%kind=1
! call kill(s2%t)
! s2%alloc=f
endif
else ! S2 does not exist
if(s2%i==0) then
S2%R=R1
s2%kind=1
elseif(s2%i>0.and.s2%i<=nv) then
call alloc(s2%t)
s2%t=(/R1,s2%S/).var.s2%i
! call var(s2%t,R1,s2%S,s2%i)
! s2%i=0
s2%kind=2
s2%alloc=t
else
line= "trouble in Dequaldacon in Real_polymorph"
ipause=mypauses(-788,line)
endif
endif ! S2 not allocated
! master=localmaster
END SUBROUTINE Dequaldacon
SUBROUTINE iequaldaconn(S2,R1)
implicit none
type (real_8),INTENT(inOUT),dimension(:)::S2
integer,INTENT(IN)::R1
integer i
! integer localmaster
! localmaster=master
! call check_snake
! 2002.12.26 master=0
do i=1,r1
call resetpoly(s2(i))
s2(i)%i=i
enddo
! master=localmaster
END SUBROUTINE iequaldaconn
SUBROUTINE equaldacon(S2,R1)
implicit none
integer ipause, mypauses
type (real_8),INTENT(inOUT)::S2
real(sp),INTENT(IN)::R1
if(real_warning) call real_stop
IF(S2%KIND==M3) THEN
if(setknob) then
s2%r=REAL(R1,kind=DP)
return
else
line= "Forbidden in equaldacon: s2 is a knob "
ipause=mypauses(0,line)
endif
ENDIF
!localmaster=master
! call check_snake
! master=0
if (s2%kind/=0) then ! S2 exist
if(S2%kind==1.or.s2%kind==3) then
S2%R=REAL(R1,kind=DP)
else
! S2%t=REAL(R1,kind=DP)
! S2%t=R1
s2%r=REAL(R1,kind=DP)
s2%kind=1
! call kill(s2%t)
! s2%alloc=f
endif
else ! S2 does not exist
if(s2%i==0) then
S2%R=REAL(R1,kind=DP)
s2%kind=1
elseif(s2%i>0.and.s2%i<=nv) then
call alloc(s2%t)
s2%t=(/REAL(R1,kind=DP),S2%S/).var.s2%i
! call var(s2%t,REAL(R1,kind=DP),S2%S,s2%i)
! s2%i=0
s2%kind=2
s2%alloc=t
else
stop 779
endif
endif ! S2 not allocated
! master=localmaster
END SUBROUTINE equaldacon
SUBROUTINE iequaldacon(S2,R1)
implicit none
integer ipause, mypauses
type (real_8),INTENT(inOUT)::S2
integer,INTENT(IN)::R1
IF(S2%KIND==M3) THEN
if(setknob) then
s2%r=r1
return
else
line= "Forbidden in iequaldacon: s2 is a knob "
ipause=mypauses(0,line)
endif
ENDIF
!localmaster=master
! call check_snake
! master=0
if (s2%kind/=0) then ! S2 exist
if(S2%kind==1.or.s2%kind==3) then
S2%R=REAL(R1,kind=DP)
else
! S2%t=REAL(R1,kind=DP)
s2%r=REAL(R1,kind=DP)
s2%kind=1
! call kill(s2%t)
! s2%alloc=f
endif
else ! S2 does not exist
if(s2%i==0) then
S2%R=REAL(R1,kind=DP)
s2%kind=1
elseif(s2%i>0.and.s2%i<=nv) then
call alloc(s2%t)
s2%t=(/REAL(R1,kind=DP),S2%S/).var.s2%i
! call var(s2%t,REAL(R1,kind=DP),S2%S,s2%i)
! s2%i=0
s2%kind=2
s2%alloc=t
else
stop 776
endif
endif ! S2 not allocated
! master=localmaster
END SUBROUTINE iequaldacon
subroutine assp(s1)
implicit none
TYPE (real_8) s1
integer ipause,mypauses
! lastmaster=master ! 2002.12.13
select case(master)
case(0:ndumt-1)
master=master+1
case(ndumt)
line= " cannot indent anymore "
ipause=mypauses(0,line)
end select
! write(26,*) " real polymorph ",master
call ass0(s1%t)
s1%alloc=t
s1%kind=2
s1%i=0
end subroutine ASSp
subroutine assp_no_master(s1)
implicit none
TYPE (real_8) s1
call ass0(s1%t)
s1%alloc=t
s1%kind=1
s1%i=0
end subroutine assp_no_master
subroutine assp_master(s1)
implicit none
TYPE (real_8) s1
integer ipause,mypauses
! lastmaster=master ! 2002.12.13
select case(master)
case(0:ndumt-1)
master=master+1
case(ndumt)
line= " cannot indent anymore "
ipause=mypauses(0,line)
end select
! write(26,*) " real polymorph ",master
call ass0(s1%t)
s1%alloc=t
s1%kind=1
s1%i=0
end subroutine assp_master
! complex
FUNCTION datant( S1 )
implicit none
TYPE (real_8) datant
TYPE (real_8), INTENT (IN) :: S1
TYPE (complextaylor) w
integer localmaster
select case(s1%kind)
case(m1)
datant%r=atan(s1%r)
datant%kind=1
case(m2)
localmaster=master
call ass(datant)
call alloc(w)
w%r=s1%t
w= atan(w)
datant%t=w%r
call kill(w)
master=localmaster
case(m3)
if(knob) then
localmaster=master
call ass(datant)
call alloc(w)
call varfk1(S1)
w%r=varf1
w= atan(w)
datant%t=w%r
call kill(w)
master=localmaster
else
datant%r=atan(s1%r)
datant%kind=1
endif
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in datant "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END FUNCTION datant
FUNCTION datanht( S1 )
implicit none
TYPE (real_8) datanht
TYPE (real_8), INTENT (IN) :: S1
integer localmaster
select case(s1%kind)
case(m1)
datanht%r=log((1+s1%r)/(1-s1%r))/2.0_dp
datanht%kind=1
case(m2)
localmaster=master
call ass(datanht)
datanht%t=atanh(s1%t)
master=localmaster
case(m3)
if(knob) then
localmaster=master
call ass(datanht)
call varfk1(S1)
datanht%t=atanh(varf1)
master=localmaster
else
datanht%r=log((1+s1%r)/sqrt(1-s1%r))/2.0_dp
datanht%kind=1
endif
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in datant "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END FUNCTION datanht
FUNCTION datanDt( S1 )
implicit none
TYPE (real_8) datanDt
TYPE (real_8), INTENT (IN) :: S1
TYPE (complextaylor) w
integer localmaster
select case(s1%kind)
case(m1)
datanDt%r=atan(s1%r)*RAD_TO_DEG_
datanDt%kind=1
case(m2)
localmaster=master
call ass(datanDt)
call alloc(w)
w%r=s1%t
w= atan(w)
datanDt%t=w%r*RAD_TO_DEG_
call kill(w)
master=localmaster
case(m3)
if(knob) then
localmaster=master
call ass(datanDt)
call alloc(w)
call varfk1(S1)
w%r=varf1
w= atan(w)
datanDt%t=w%r*RAD_TO_DEG_
call kill(w)
master=localmaster
else
datanDt%r=atan(s1%r)*RAD_TO_DEG_
datanDt%kind=1
endif
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in datanDt "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END FUNCTION datanDt
FUNCTION absoftdatanDr( S1 )
implicit none
real(dp) absoftdatanDr
real(dp), INTENT (IN) :: S1
absoftdatanDr=atan(s1)*RAD_TO_DEG_
END FUNCTION absoftdatanDr
FUNCTION datan2t( S2,S1 )
implicit none
integer ipause, mypauses
TYPE (real_8) datan2t
TYPE (real_8), INTENT (IN) :: S2,S1
TYPE (complextaylor) w
real(dp) ANG
integer localmaster ,si
si=1.0_dp
if(s2%kind==0.or.S1%kind==0) then
line= " Problems in datan2t "
ipause=mypauses(0,line)
endif
if(S2%kind==S1%kind) then ! 1
select case(S1%kind)
case(m1)
DATAN2T%R=ATAN2(S2%R,S1%R)
datan2t%kind=1
case(m2)
localmaster=master
call ass(datan2t)
call alloc(w)
ANG=ATAN2(S2%T.SUB.'0',S1%T.SUB.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then
w%r=s2%t/s1%t
w=atan(w)
datan2t%t=w%r
datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
else
if((s2%t.sub.'0')<0.0_dp) si=-1.0_dp
w%r=s1%t/SQRT(s2%t**2+s1%t**2)
w=acos(w)
datan2t%t=w%r
datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
endif
call kill(w)
master=localmaster
case(m3)
localmaster=master
if(knob) then !knob 1
call ass(datan2t)
call alloc(w)
call varfk1(S1)
call varfk2(S2)
ANG=ATAN2(varf2.SUB.'0',varf1.SUB.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then
w%r=varf2/varf1
w=atan(w)
datan2t%t=w%r
datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
else
if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
w%r=varf1/SQRT(varf2**2+varf1**2)
w=acos(w)
datan2t%t=w%r
datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
endif
call kill(w)
else !knob 1
DATAN2T%R=ATAN2(S2%R,S1%R)
datan2t%kind=1
endif !knob 1
master=localmaster
end select
elseif(S2%kind>S1%kind) then !1
localmaster=master
! call ass(datan2t)
call alloc(w)
if(s2%kind==2) then !
call ass(datan2t)
ANG=ATAN2(S2%T.SUB.'0',S1%R)
if(pil>abs(ANG).or.abs(ANG)>pim) then
w%r=s2%t/s1%r
w=atan(w)
datan2t%t=w%r
datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
else
if((s2%t.sub.'0')<0.0_dp) si=-1.0_dp
w%r=s1%t/SQRT(s2%t**2+s1%r**2)
w=acos(w)
datan2t%t=w%r
datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
endif
else !
if(knob) then ! knob
call ass(datan2t)
call varfk2(S2)
if(s1%kind==1) then !!
ANG=ATAN2(varf2.SUB.'0',S1%R)
if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
w%r=varf2/s1%r
w=atan(w)
datan2t%t=w%r
datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
else !!!
if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
w%r=s1%R/SQRT(varf2**2+s1%R**2) ! etienne error here????
w=acos(w)
datan2t%t=w%r
datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
endif !!!
else !!
ANG=ATAN2(varf2.SUB.'0',S1%t.sub.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
w%r=varf2/s1%t
w=atan(w)
datan2t%t=w%r
datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
else !!!
if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
w%r=s1%t/SQRT(varf2**2+s1%t**2)
w=acos(w)
datan2t%t=w%r
datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
endif !!!
endif !!
else ! knob
if(s1%kind==1) then !!
DATAN2T%R=ATAN2(S2%R,S1%R)
datan2t%kind=1
else !!
call ass(datan2t)
ANG=ATAN2(S2%R,S1%t.sub.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
w%r=S2%R/s1%t
w=atan(w)
datan2t%t=w%r
datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
else !!!
if((S2%R)<0.0_dp) si=-1.0_dp
w%r=s1%t/SQRT(S2%R**2+s1%t**2)
w=acos(w)
datan2t%t=w%r
datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
endif !!!
endif !!
endif ! knob
endif !
call kill(w)
master=localmaster
ELSE ! 1
localmaster=master
! call ass(datan2t)
call alloc(w)
if(s1%kind==2) then !
call ass(datan2t)
ANG=ATAN2(S2%R,S1%T.SUB.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then
w%r=s2%r/s1%t
w=atan(w)
datan2t%t=w%r
datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
else
if(s2%r<0.0_dp) si=-1.0_dp
w%r=s1%t/SQRT(s2%r**2+s1%t**2)
w=acos(w)
datan2t%t=w%r
datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
endif
else !
IF(KNOB) THEN ! KNOB 2
call ass(datan2t)
if(s2%kind==1) then !!
call varfk1(S1)
ANG=ATAN2(S2%R,varf1.SUB.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
w%r=s2%r/varf1
w=atan(w)
datan2t%t=w%r
datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
else !!!
if(s2%r<0.0_dp) si=-1.0_dp
w%r=varf1/SQRT(s2%r**2+varf1**2)
w=acos(w)
datan2t%t=w%r
datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
endif !!!
else !! s2%kind ==2
call varfk1(S1)
ANG=ATAN2(S2%t.sub.'0',varf1.SUB.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
w%r=s2%t/varf1
w=atan(w)
datan2t%t=w%r
datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
else !!!
if((S2%t.sub.'0')<0.0_dp) si=-1.0_dp
w%r=varf1/SQRT(s2%t**2+varf1**2)
w=acos(w)
datan2t%t=w%r
datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
endif !!!
endif !!
ELSE !KNOB 2
if(s2%kind==1) then !!
DATAN2T%R=ATAN2(S2%R,S1%R)
datan2t%kind=1
else !! s2%kind ==2
call ass(datan2t)
ANG=ATAN2(S2%t.sub.'0',S1%r)
if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
w%r=s2%t/S1%r
w=atan(w)
datan2t%t=w%r
datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
else !!!
if((S2%t.sub.'0')<0.0_dp) si=-1.0_dp
w%r=S1%r/SQRT(s2%t**2+S1%r**2)
w=acos(w)
datan2t%t=w%r
datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
endif !!!
endif !!
ENDIF ! KNOB 2
endif !
call kill(w)
master=localmaster
endif
END FUNCTION datan2t
FUNCTION datan2dt( S2,S1 )
implicit none
integer ipause, mypauses
TYPE (real_8) datan2dt
TYPE (real_8), INTENT (IN) :: S2,S1
TYPE (complextaylor) w
real(dp) ANG
integer localmaster ,si
si=1.0_dp
if(s2%kind==0.or.S1%kind==0) then
line= " Problems in datand2t "
ipause=mypauses(0,line)
endif
if(S2%kind==S1%kind) then
select case(S1%kind)
case(m1)
datan2dt%R=ATAN2(S2%R,S1%R)*RAD_TO_DEG_
datan2dt%kind=1
case(m2)
localmaster=master
call ass(datan2dt)
call alloc(w)
ANG=ATAN2(S2%T.SUB.'0',S1%T.SUB.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then
w%r=s2%t/s1%t
w=atan(w)
datan2dt%t=w%r
datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
else
if((s2%t.sub.'0')<0.0_dp) si=-1.0_dp
w%r=s1%t/SQRT(s2%t**2+s1%t**2)
w=acos(w)
datan2dt%t=w%r
datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
endif
call kill(w)
master=localmaster
case(m3)
localmaster=master
if(knob) then ! knob 1
call ass(datan2dt)
call alloc(w)
call varfk1(S1)
call varfk2(S2)
ANG=ATAN2(varf2.SUB.'0',varf1.SUB.'0') !etienne
if(pil>abs(ANG).or.abs(ANG)>pim) then
w%r=varf2/varf1
w=atan(w)
datan2dt%t=w%r
datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
else
if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
w%r=varf1/SQRT(varf2**2+varf1**2)
w=acos(w)
datan2dt%t=w%r
datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
endif
call kill(w)
else ! knob 1
datan2dt%R=ATAN2(S2%R,S1%R)*RAD_TO_DEG_
datan2dt%kind=1
endif ! knob 1
master=localmaster
end select
elseif(S2%kind>S1%kind) then
localmaster=master
! call ass(datan2dt)
call alloc(w)
if(s2%kind==2) then !
call ass(datan2dt)
ANG=ATAN2(S2%T.SUB.'0',S1%R)
if(pil>abs(ANG).or.abs(ANG)>pim) then
w%r=s2%t/s1%r
w=atan(w)
datan2dt%t=w%r
datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
else
if((s2%t.sub.'0')<0.0_dp) si=-1.0_dp
w%r=s1%t/SQRT(s2%t**2+s1%r**2)
w=acos(w)
datan2dt%t=w%r
datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
endif
else !
if(knob) then !knob 2
call ass(datan2dt)
call varfk2(S2)
if(s1%kind==1) then !!
ANG=ATAN2(varf2.SUB.'0',S1%R)
if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
w%r=varf2/s1%r
w=atan(w)
datan2dt%t=w%r
datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
else !!!
if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
w%r=s1%t/SQRT(varf2**2+s1%t**2)
w=acos(w)
datan2dt%t=w%r
datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
endif !!!
else !!
ANG=ATAN2(varf2.SUB.'0',S1%t.sub.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
w%r=varf2/s1%t
w=atan(w)
datan2dt%t=w%r
datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
else !!!
if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
w%r=s1%t/SQRT(varf2**2+s1%t**2)
w=acos(w)
datan2dt%t=w%r
datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
endif !!!
endif !!
else !knob 2
if(s1%kind==1) then !!
datan2dt%R=ATAN2(S2%R,S1%R)*RAD_TO_DEG_
datan2dt%kind=1
else !!
call ass(datan2dt)
ANG=ATAN2(S2%R,S1%t.sub.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
w%r=S2%R/s1%t
w=atan(w)
datan2dt%t=w%r
datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
else !!!
if((S2%R)<0.0_dp) si=-1.0_dp
w%r=s1%t/SQRT(S2%R**2+s1%t**2)
w=acos(w)
datan2dt%t=w%r
datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
endif !!!
endif !!
endif !knob 2
endif !
call kill(w)
master=localmaster
ELSE !1
localmaster=master
! call ass(datan2dt)
call alloc(w)
if(s1%kind==2) then !
call ass(datan2dt)
ANG=ATAN2(S2%R,S1%T.SUB.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then
w%r=s2%r/s1%t
w=atan(w)
datan2dt%t=w%r
datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
else
if(s2%r<0.0_dp) si=-1.0_dp
w%r=s1%t/SQRT(s2%r**2+s1%t**2)
w=acos(w)
datan2dt%t=w%r
datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
endif
else !
if(knob) then ! knob 3
call ass(datan2dt)
if(s2%kind==1) then !!
call varfk1(S1)
ANG=ATAN2(S2%R,varf1.SUB.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then
w%r=s2%r/varf1
w=atan(w)
datan2dt%t=w%r
datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
else
if(s2%r<0.0_dp) si=-1.0_dp
w%r=varf1/SQRT(s2%r**2+varf1**2)
w=acos(w)
datan2dt%t=w%r
datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
endif
else !! s2%kind ==2
call varfk1(S1)
ANG=ATAN2(S2%t.sub.'0',varf1.SUB.'0')
if(pil>abs(ANG).or.abs(ANG)>pim) then
w%r=s2%T/varf1
w=atan(w)
datan2dt%t=w%r
datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
else
if((S2%t.sub.'0')<0.0_dp) si=-1.0_dp
w%r=varf1/SQRT(s2%T**2+varf1**2)
w=acos(w)
datan2dt%t=w%r
datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
endif
endif !!
else ! knob 3
if(s2%kind==1) then !!
datan2dt%R=ATAN2(S2%R,S1%R)*RAD_TO_DEG_
datan2dt%kind=1
else !! s2%kind ==2
call ass(datan2dt)
ANG=ATAN2(S2%t.sub.'0',S1%r)
if(pil>abs(ANG).or.abs(ANG)>pim) then
w%r=s2%T/S1%r
w=atan(w)
datan2dt%t=w%r
datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
else
if((S2%t.sub.'0')<0.0_dp) si=-1.0_dp
w%r=S1%r/SQRT(s2%T**2+S1%r**2)
w=acos(w)
datan2dt%t=w%r
datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
datan2dt%t=datan2dt%t*RAD_TO_DEG_
endif
endif !!
endif ! knob 3
endif !
call kill(w)
master=localmaster
endif
END FUNCTION datan2dt
FUNCTION absoftdatan2dr( S2,S1 )
implicit none
real(dp) absoftdatan2dr
real(dp), INTENT (IN) :: S2,S1
absoftdatan2dr=ATAN2(S2,S1)*RAD_TO_DEG_
END FUNCTION absoftdatan2dr
FUNCTION dasint( S1 )
implicit none
TYPE (real_8) dasint
TYPE (real_8), INTENT (IN) :: S1
TYPE (taylor) w
integer localmaster
select case(s1%kind)
case(m1)
dasint%r=arcsin(s1%r)
! dasint%r=asin(s1%r)
dasint%kind=1
case(m2)
localmaster=master
call ass(dasint)
call alloc(w)
w=s1%t
w= asin(w)
dasint%t=w
call kill(w)
master=localmaster
case(m3)
if(knob) then
localmaster=master
call ass(dasint)
call alloc(w)
call varfk1(S1)
w=varf1
w= asin(w)
dasint%t=w
call kill(w)
master=localmaster
else
dasint%r=asin(s1%r)
dasint%kind=1
endif
! call unvarkind3(S1)
end select
END FUNCTION dasint
FUNCTION dacost( S1 )
implicit none
TYPE (real_8) dacost
TYPE (real_8), INTENT (IN) :: S1
TYPE (taylor) w
integer localmaster
select case(s1%kind)
case(m1)
dacost%r=arccos(s1%r)
dacost%kind=1
case(m2)
localmaster=master
call ass(dacost)
call alloc(w)
w=s1%t
w= acos(w)
dacost%t=w
call kill(w)
master=localmaster
case(m3)
if(knob) then
localmaster=master
call ass(dacost)
call alloc(w)
call varfk1(S1)
w=varf1
w= acos(w)
dacost%t=w
call kill(w)
master=localmaster
else
dacost%r=acos(s1%r)
dacost%kind=1
endif
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in dacost "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END FUNCTION dacost
FUNCTION dcosht( S1 )
implicit none
TYPE (real_8) dcosht
TYPE (real_8), INTENT (IN) :: S1
TYPE (complextaylor) w
integer localmaster
select case(s1%kind)
case(m1)
dcosht%r=cosh(s1%r)
dcosht%kind=1
case(m2)
localmaster=master
call ass(dcosht)
call alloc(w)
w%r=s1%t
w= cosh(w)
dcosht%t=w%r
call kill(w)
master=localmaster
case(m3)
if(knob) then
localmaster=master
call ass(dcosht)
call alloc(w)
call varfk1(S1)
w%r=varf1
w= cosh(w)
dcosht%t=w%r
call kill(w)
master=localmaster
else
dcosht%r=cosh(s1%r)
dcosht%kind=1
endif
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in dcosht "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END FUNCTION dcosht
FUNCTION dsinht( S1 )
implicit none
TYPE (real_8) dsinht
TYPE (real_8), INTENT (IN) :: S1
TYPE (complextaylor) w
integer localmaster
select case(s1%kind)
case(m1)
dsinht%r=sinh(s1%r)
dsinht%kind=1
case(m2)
localmaster=master
call ass(dsinht)
call alloc(w)
w%r=s1%t
w= sinh(w)
dsinht%t=w%r
call kill(w)
master=localmaster
case(m3)
if(knob) then
localmaster=master
call ass(dsinht)
call alloc(w)
call varfk1(S1)
w%r=varf1
w= sinh(w)
dsinht%t=w%r
call kill(w)
master=localmaster
else
dsinht%r=sinh(s1%r)
dsinht%kind=1
endif
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in dsinht "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END FUNCTION dsinht
FUNCTION dtanht( S1 )
implicit none
TYPE (real_8) dtanht
TYPE (real_8), INTENT (IN) :: S1
integer localmaster
select case(s1%kind)
case(m1)
dtanht%r=tanh(s1%r)
dtanht%kind=1
case(m2)
localmaster=master
call ass(dtanht)
dtanht%t=tanh(s1%t)
master=localmaster
case(m3)
if(knob) then
localmaster=master
call ass(dtanht)
call varfk1(S1)
dtanht%t=tanh(varf1)
master=localmaster
else
dtanht%r=tanh(s1%r)
dtanht%kind=1
endif
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in dtanht "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END FUNCTION dtanht
SUBROUTINE MAPreal_8(S2,S1)
implicit none
type (damap),INTENT(inout)::S2
type (real_8),INTENT(IN),dimension(:)::S1
integer i
call check_snake
do i=1,nd2
s2%v(i)=s1(i) !%t
enddo
END SUBROUTINE MAPreal_8
SUBROUTINE normal_p(S2,S1)
implicit none
type (normalform),INTENT(inOUT)::S2
type (real_8),INTENT(IN),dimension(:)::S1
type (damap) t
integer i
call alloc(t)
call check_snake
do i=1,nd2
t%v(i)=s1(i) !%t
enddo
s2=t
call kill(t)
END SUBROUTINE normal_p
SUBROUTINE real_8MAP(S1,S2)
implicit none
type (damap),INTENT(in)::S2
type (real_8),INTENT(inout),dimension(:)::S1
integer i
call check_snake
do i=1,nd2
s1(i)=s2%v(i)
enddo
END SUBROUTINE real_8MAP
! radiation
SUBROUTINE ENVreal_8(S2,S1)
implicit none
type (ENV_8),INTENT(inout),dimension(:)::S2
type (real_8),INTENT(IN),dimension(:)::S1
integer i
call check_snake
do i=1,nd2
s2(i)%v=s1(i) !%t
enddo
END SUBROUTINE ENVreal_8
SUBROUTINE real_8ENV(S1,S2)
implicit none
type (ENV_8),INTENT(in),dimension(:)::S2
type (real_8),INTENT(inout),dimension(:)::S1
integer i
call check_snake
do i=1,nd2
s1(i)=s2(i)%v
enddo
END SUBROUTINE real_8ENV
SUBROUTINE ENV_8RADTAYLOR(S2,S1)
implicit none
type (ENV_8),INTENT(inout),dimension(:)::S2
type (RADTAYLOR),INTENT(IN),dimension(:)::S1
integer i,J
call check_snake
do i=1,nd2
s2(i)%v=s1(i)%V !%t
enddo
do i=1,nd2
do J=1,nd2
s2(i)%E(J)=s1(i)%E(J) !%t
enddo
enddo
END SUBROUTINE ENV_8RADTAYLOR
SUBROUTINE RADTAYLORENV_8(S1,S2)
implicit none
type (ENV_8),INTENT(in),dimension(:)::S2
type (RADTAYLOR),INTENT(inout),dimension(:)::S1
integer i,J
call check_snake
do i=1,nd2
s1(i)%V= s2(i)%v !%t
enddo
do i=1,nd2
do J=1,nd2
s1(i)%E(J)=s2(i)%E(J) !%t
enddo
enddo
END SUBROUTINE RADTAYLORENV_8
SUBROUTINE RADTAYLORprobe_8(S1,S2)
implicit none
type (probe_8),INTENT(in) ::S2
type (RADTAYLOR),INTENT(inout),dimension(:)::S1
integer i,J
call check_snake
do i=1,nd2
s1(i)%V= s2%X(I) !%t
enddo
do i=1,nd2
do J=1,nd2
s1(i)%E(J)=s2%E_IJ(i,J) !%t
enddo
enddo
END SUBROUTINE RADTAYLORprobe_8
SUBROUTINE beamprobe_8(S2,Sprobe_8)
implicit none
type (beamenvelope),INTENT(inout)::S2
type (radtaylor) S1(ndim2)
type (probe_8),INTENT(IN)::Sprobe_8
call alloc(s1,nd2)
call check_snake
S1=Sprobe_8
s2=s1
call kill(s1,nd2)
END SUBROUTINE beamprobe_8
SUBROUTINE resetenv(S2)
implicit none
type (env_8),INTENT(INOUT)::S2
integer i
call resetpoly(s2%v)
do i=1,nd2
call resetpoly(s2%e(i))
call resetpoly(s2%sigma0(i))
call resetpoly(s2%sigmaf(i))
enddo
END SUBROUTINE resetenv
SUBROUTINE allocenv(S2)
implicit none
type (env_8),INTENT(INOUT)::S2
integer i
call alloc(s2%v)
do i=1,nd2
call alloc(s2%e(i))
call alloc(s2%sigma0(i))
call alloc(s2%sigmaf(i))
enddo
END SUBROUTINE allocenv
SUBROUTINE allocenvn(S2,K)
implicit none
type (env_8),INTENT(INOUT),dimension(:)::S2
INTEGER,optional,INTENT(IN)::k
INTEGER J,i,N
if(present(k)) then
I=LBOUND(S2,DIM=1)
N=LBOUND(S2,DIM=1)+K-1
else
I=LBOUND(S2,DIM=1)
N=UBOUND(S2,DIM=1)
endif
DO J=I,N
call allocenv(s2(j))
enddo
END SUBROUTINE allocenvn
SUBROUTINE e_OPT(S1,S2,s3,s4,s5,s6,s7,s8,s9,s10)
implicit none
type (env_8),INTENT(INout)::S1
type (env_8),optional, INTENT(INout):: S2,s3,s4,s5,s6,s7,s8,s9,s10
call allocenv(s1)
if(present(s2)) call allocenv(s2)
if(present(s3)) call allocenv(s3)
if(present(s4)) call allocenv(s4)
if(present(s5)) call allocenv(s5)
if(present(s6)) call allocenv(s6)
if(present(s7)) call allocenv(s7)
if(present(s8)) call allocenv(s8)
if(present(s9)) call allocenv(s9)
if(present(s10))call allocenv(s10)
END SUBROUTINE e_opt
! SUBROUTINE allocenvn(S2,n)
! implicit none
! type (env_8),INTENT(INOUT),dimension(:)::S2
! integer, INTENT(IN)::n
! integer i
!
! do i=1,n
! call allocenv(s2(i))
! enddo
! END SUBROUTINE allocenvn
SUBROUTINE resetenvn(S2,K)
implicit none
type (env_8),INTENT(INOUT),dimension(:)::S2
INTEGER,optional,INTENT(IN)::k
INTEGER J,i,N
if(present(k)) then
I=LBOUND(S2,DIM=1)
N=LBOUND(S2,DIM=1)+K-1
else
I=LBOUND(S2,DIM=1)
N=UBOUND(S2,DIM=1)
endif
DO J=I,N
call resetenv(s2(j))
enddo
END SUBROUTINE resetenvn
! SUBROUTINE resetenvn(S2,n)
! implicit none
! type (env_8),INTENT(INOUT),dimension(:)::S2
! integer, INTENT(IN)::n
! integer i
!
! do i=1,n
! call resetenv(s2(i))
! enddo
!
!
! END SUBROUTINE resetenvn
SUBROUTINE beamENV_8(S2,SENV_8)
implicit none
type (beamenvelope),INTENT(inout)::S2
type (radtaylor) S1(ndim2)
! type (ENV_8),INTENT(IN)::SENV_8(ndim2)
type (ENV_8),INTENT(IN)::SENV_8(6)
call alloc(s1,nd2)
call check_snake
S1=SENV_8
s2=s1
call kill(s1,nd2)
END SUBROUTINE beamENV_8
SUBROUTINE varfk1(S2)
implicit none
! type (real_8),INTENT(INOUT)::S2
type (real_8),INTENT(IN):: S2
if(knob) then
if(nb_==0) then
varf1=(/S2%R,S2%S/).var.(s2%i+npara_fpp)
elseif(s2%nb==nb_) then
varf1=(/S2%R,S2%S/).var.(s2%i+npara_fpp-s2%g+1)
else
varf1=S2%R
endif
else ! Not a knob
stop 333
varf1=(/S2%R,0.0_dp/).var.0 ! this is a buggy line never used
endif
end SUBROUTINE varfk1
SUBROUTINE varfk2(S2)
implicit none
! type (real_8),INTENT(INOUT)::S2
type (real_8),INTENT(IN):: S2
if(knob) then
if(nb_==0) then
varf2=(/S2%R,S2%S/).var.(s2%i+npara_fpp)
elseif(s2%nb==nb_) then
varf2=(/S2%R,S2%S/).var.(s2%i+npara_fpp-s2%g+1)
else
varf2=S2%R
endif
else ! Not a knob
stop 334
varf2=(/S2%R,0.0_dp/).var.0 ! this is a buggy line never used
endif
end SUBROUTINE varfk2
! SINH(X)/X DONE partly COSY-INFINITY WISE FOR TPSA
function SINH_HR(X)
implicit none
real(dp), INTENT (IN) :: X
real(dp) SINH_HR,Y,NORM0,NORM,SINH_HR0
logical(lp) NOTDONE,CHECK
INTEGER I,ipause,mypauses
if(abs(x)=NORM0) THEN
NOTDONE=.FALSE.
ELSE
NORM0=NORM
ENDIF
ENDIF
SINH_HR=SINH_HR0
I=I+2
ENDDO
IF(I==NMAX_pol) THEN
line="NO CONVERGENCE IN SINH_HR"
ipause=mypauses(NMAX_pol,line)
ENDIF
else
SINH_HR=SINH(X)/X
endif
return
end function SINH_HR
! SIN(X)/X DONE partly COSY-INFINITY WISE FOR TPSA
function SIN_HR(X)
implicit none
real(dp), INTENT (IN) :: X
real(dp) SIN_HR,Y,NORM0,NORM,SINH_HR0
logical(lp) NOTDONE,CHECK
INTEGER I,ipause,mypauses
if(abs(x)=NORM0) THEN
NOTDONE=.FALSE.
ELSE
NORM0=NORM
ENDIF
ENDIF
SIN_HR=SINH_HR0
I=I+2
ENDDO
IF(I==NMAX_pol) THEN
line="NO CONVERGENCE IN SINH_HR"
ipause=mypauses(NMAX_pol,line)
ENDIF
else
SIN_HR=SIN(X)/X
endif
return
end function SIN_HR
FUNCTION sinX_Xt( S1 )
implicit none
integer ipause, mypauses
TYPE (real_8) sinX_Xt
TYPE (real_8), INTENT (IN) :: S1
integer localmaster
TYPE(TAYLOR)SIN_HP
TYPE(TAYLOR)Y,SINH_HP0
real(dp) NORM0,NORM
logical(lp) NOTDONE,CHECK
INTEGER I
select case(s1%kind)
case(m1)
sinX_Xt%r=SIN_HR(s1%r)
sinX_Xt%kind=1
case(m2)
CALL ALLOC(SIN_HP); CALL ALLOC(Y); CALL ALLOC(SINH_HP0);
localmaster=master
call ass(sinX_Xt)
if(ABS(S1%T.SUB.'0')=NORM0.and.(.not.check)) THEN ! sagan
NOTDONE=.FALSE.
ELSE
NORM0=NORM
ENDIF
ENDIF
SIN_HP=SINH_HP0
I=I+2
ENDDO
IF(I==NMAX_pol) THEN
line="NO CONVERGENCE IN SIN_HP"
ipause=mypauses(NMAX_pol,line)
ENDIF
else
SIN_HP=SIN(S1%T)/S1%T
endif
sinX_Xt%t= SIN_HP
master=localmaster
CALL KILL(SIN_HP); CALL KILL(Y); CALL KILL(SINH_HP0);
case(m3)
if(knob) then
localmaster=master
call ass(sinX_Xt)
call varfk1(S1)
if(ABS(varf1.SUB.'0')=NORM0.and.(.not.check)) THEN ! sagan
NOTDONE=.FALSE.
ELSE
NORM0=NORM
ENDIF
ENDIF
SIN_HP=SINH_HP0
I=I+2
ENDDO
IF(I==NMAX_pol) THEN
line="NO CONVERGENCE IN SINH_HP"
ipause=mypauses(NMAX_pol,line)
ENDIF
else
SIN_HP=SIN(varf1)/varf1
endif
sinX_Xt%t= SIN_HP
master=localmaster
else
sinX_Xt%r= SIN_HR(S1%r)
sinX_Xt%kind=1
endif
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in sinX_XT "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END FUNCTION sinX_XT
FUNCTION sinHX_Xt( S1 )
implicit none
integer ipause, mypauses
TYPE (real_8) sinHX_Xt
TYPE (real_8), INTENT (IN) :: S1
integer localmaster
TYPE(TAYLOR)SIN_HP
TYPE(TAYLOR)Y,SINH_HP0
real(dp) NORM0,NORM
logical(lp) NOTDONE,CHECK
INTEGER I
select case(s1%kind)
case(m1)
sinHX_Xt%r=SINH_HR(s1%r)
sinHX_Xt%kind=1
case(m2)
CALL ALLOC(SIN_HP); CALL ALLOC(Y); CALL ALLOC(SINH_HP0);
localmaster=master
call ass(sinHX_Xt)
if(ABS(S1%T.SUB.'0')=NORM0.and.(.not.check)) THEN ! sagan
NOTDONE=.FALSE.
ELSE
NORM0=NORM
ENDIF
ENDIF
SIN_HP=SINH_HP0
I=I+2
ENDDO
IF(I==NMAX_pol) THEN
line="NO CONVERGENCE IN SINH_HP"
ipause=mypauses(NMAX_pol,line)
ENDIF
else
SIN_HP=SINH(S1%T)/S1%T
endif
sinHX_Xt%t= SIN_HP
master=localmaster
CALL KILL(SIN_HP); CALL KILL(Y); CALL KILL(SINH_HP0);
case(m3)
if(knob) then
localmaster=master
call ass(sinHX_Xt)
call varfk1(S1)
if(ABS(varf1.SUB.'0')=NORM0.and.(.not.check)) THEN ! sagan
NOTDONE=.FALSE.
ELSE
NORM0=NORM
ENDIF
ENDIF
SIN_HP=SINH_HP0
I=I+2
ENDDO
IF(I==NMAX_pol) THEN
line="NO CONVERGENCE IN SINH_HP"
ipause=mypauses(NMAX_pol,line)
ENDIF
else
SIN_HP=SINH(varf1)/varf1
endif
sinHX_Xt%t= SIN_HP
master=localmaster
else
sinHX_Xt%r= SINH_HR(S1%r)
sinHX_Xt%kind=1
endif
case default
w_p=0
w_p%nc=2
w_p%fc='((1X,A72,/,1x,a72))'
w_p%fi='(2((1X,i4)))'
w_p%c(1)= " trouble in sinHX_Xt "
w_p%c(2)= "s1%kind "
w_p=(/s1%kind /)
! call !write_e(0)
end select
END FUNCTION sinHX_Xt
! remove small numbers
SUBROUTINE clean_real_8(S1,S2,prec)
implicit none
type (real_8),INTENT(INOUT)::S2
type (real_8), intent(INOUT):: s1
real(dp) prec
integer i
type(real_8) t
call alloc(t)
t=s1
select case(s1%kind)
case(m1)
if(abs(t%r)