C### DAEPACK v0.02 - Copyright (C) M.I.T. C### - DISCONTINUITY LOCKED MODEL - subroutine resdl(neq,t,w,wdot,delta,ires,ichvar,rpar,ipar, $ zzznsc,zzzndf,zzzisc,zzziscchng,zzzdiscon,zzzlocked, $ zzzistate,zzzig) C### implicit none integer neq,ires,ichvar,ipar(1) double precision t,w(neq),wdot(neq),delta(neq),rpar(1) double precision x1,x2,x3,x1dot,x2dot,v1,v2,v3,z3,z4 double precision a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,a11,a12 parameter(a1=12.0d0,a2=0.24d0,a3=12.0d0,a4=0.24d0) parameter(a5=13.64d0,a6=-11.64d0,a7=-50.0d0,a8=0.0d0) parameter(a9=-11.64d0,a10=13.64d0,a11=0.0d0,a12=-50.0d0) double precision r3,l3 parameter(r3=10.0d0,l3=0.2d0) double precision cos,atan intrinsic cos,atan C### C### Additional argument for discontinuity-locked evaluations C### zzznsc - Total number of currently active state conditions. C### zzzndf - Total number of currently active discontinuity C### functions. C### zzzisc(zzznsc) - Number of discontinuity functions in each state C### condition. C### zzziscchng(zzznsc) - Contains indices of state conditions that have C### flipped. C### zzzdiscon(zzzndf) - Current values of discontinuity functions. C### zzzlocked - Flag indicating whether or not this is a locked C### evaluation. C### zzzistate - State index to determine whether or not it has C### flipped. C### zzzig(zzzndf) - Values of discontinuity functions to determine C### whether or not state condition has flipped. integer zzznsc,zzzndf,zzzisc(zzznsc),zzziscchng(zzznsc), integer zzzistate,zzzig(zzzndf),izzzchng double precision zzzdiscon(zzzndf) logical zzzlogval,zzzlocked,lzzzdlk external lzzzdlk,izzzchng C### zzznsc=0 zzzndf=0 C### x1dot=wdot(1) x2dot=wdot(2) x1=w(1) x2=w(2) x3=w(3) v1=w(4) v2=w(5) v3=w(6) z3=w(7) z4=w(8) delta(1)=x1+x2-x3 delta(2)=r3*(x1+x2)+l3*(x1dot+x2dot)-v3 delta(3)=100.0d0*cos(100.0d0*atan(1.0d0)*4.0d0*t)-v1 delta(4)=-v1-v2 delta(5)=v1-v3-z3 delta(6)=v2-v3-z4 C### Modified if equation for discontinuity locking zzznsc=zzznsc+1 if(zzznsc.eq.zzzistate) then zzzlogval=(zzzig(1).eq.1.or.zzzig(2).eq.1).and. $ (zzzig(3).eq.1.and.zzzig(4).eq.1) zzzistate=izzzchng(zzzistate,zzzlogval) return else zzzndf=zzzndf+1 zzzzzzdiscon(zzzndf)=x1-(0.0d0) zzzndf=zzzndf+1 zzzzzzdiscon(zzzndf)=v1-(v3) zzzndf=zzzndf+1 zzzzzzdiscon(zzzndf)=0.0d0-(x2) zzzndf=zzzndf+1 zzzzzzdiscon(zzzndf)=v3-(v2) if(.not.zzzlocked) zzzisc(zzznsc)=4 end if zzzlogval=(x1.gt.0.0d0.or.v1.gt.v3).and. $ (x2.le.0.0d0.and.v2.lt.v3) if(lzzzdlk(zzzlocked,zzznsc,zzzlogval,zzziscchng)) then delta(7)=(v1-a1*x1)/a2-x1dot delta(8)=0.0d0-x2dot else C### Modified if equation for discontinuity locking zzznsc=zzznsc+1 if(zzznsc.eq.zzzistate) then zzzlogval=(zzzig(1).eq.1.or.zzzig(2).eq.1).and. $ (zzzig(3).eq.1.and.zzzig(4).eq.1) zzzistate=izzzchng(zzzistate,zzzlogval) return else zzzndf=zzzndf+1 zzzzzzdiscon(zzzndf)=x2-(0.0d0) zzzndf=zzzndf+1 zzzzzzdiscon(zzzndf)=v2-(v3) zzzndf=zzzndf+1 zzzzzzdiscon(zzzndf)=0.0d0-(x1) zzzndf=zzzndf+1 zzzzzzdiscon(zzzndf)=v3-(v1) if(.not.zzzlocked) zzzisc(zzznsc)=4 end if zzzlogval=(x2.gt.0.0d0.or.v2.gt.v3).and. $ (x1.lt.0.0d0.and.v1.lt.v3) if(lzzzdlk(zzzlocked,zzznsc,zzzlogval,zzziscchng)) then delta(7)=0.0d0-x1dot delta(8)=(v2-a3*x2)/a4-x2dot else delta(7)=a5*v1+a6*v2+a7*x1+a8*x2-x1dot delta(8)=a9*v1+a10*v2+a11*x1+a12*x2-x2dot end if end if return C### Done end C###