File: cir.f90

package info (click to toggle)
auto-07p 0.9.1%2Bdfsg-7
  • links: PTS, VCS
  • area: main
  • in suites: buster
  • size: 16,200 kB
  • sloc: fortran: 22,644; f90: 19,340; python: 19,045; ansic: 11,116; sh: 1,079; makefile: 618; perl: 339
file content (145 lines) | stat: -rw-r--r-- 3,487 bytes parent folder | download | duplicates (5)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
!---------------------------------------------------------------------- 
!---------------------------------------------------------------------- 
!   cir :    Homoclinic Bifurcation in an Electronic Circuit
!                (the same equations as in demo tor)
!---------------------------------------------------------------------- 
!---------------------------------------------------------------------- 

      SUBROUTINE FUNC(NDIM,U,ICP,PAR,IJAC,F,DFDU,DFDP) 
!     ---------- ---- 

      IMPLICIT NONE
      INTEGER, INTENT(IN) :: NDIM, ICP(*), IJAC
      DOUBLE PRECISION, INTENT(IN) :: U(NDIM), PAR(*)
      DOUBLE PRECISION, INTENT(OUT) :: F(NDIM)
      DOUBLE PRECISION, INTENT(INOUT) :: DFDU(NDIM,NDIM), DFDP(NDIM,*)

      DOUBLE PRECISION rn,be,ga,r,a3,b3,x,y,z

       rn=PAR(1) 
       be=PAR(2) 
       ga=PAR(3) 
       r =PAR(4) 
       a3=PAR(5) 
       b3=PAR(6) 

       x=U(1)
       y=U(2)
       z=U(3)

       F(1)= ( -(be+rn)*x + be*y - a3*x**3 + b3*(y-x)**3 )/r
       F(2)=  be*x - (be+ga)*y - z - b3*(y-x)**3
       F(3)= y

       IF(IJAC.EQ.0)RETURN

       DFDU(1,1)=( -(be+rn) -3*a3*x**2 - 3*b3*(y-x)**2  )/r
       DFDU(1,2)=( be + 3*b3*(y-x)**2 )/r
       DFDU(1,3)=0

       DFDU(2,1)=be + 3*b3*(y-x)**2
       DFDU(2,2)=-(be+ga) - 3*b3*(y-x)**2
       DFDU(2,3)=-1

       DFDU(3,1)=0
       DFDU(3,2)=1
       DFDU(3,3)=0

      IF(IJAC.EQ.1)RETURN 

!      *Parameter derivatives
       DFDP(1,1)=-x/r
       DFDP(2,1)=0
       DFDP(3,1)=0

       DFDP(1,2)=( -x + y )/r
       DFDP(2,2)=x-y
       DFDP(3,2)=0

       DFDP(1,3)=0
       DFDP(2,3)=-y
       DFDP(3,3)=0

       DFDP(1,4)=-F(1)/r
       DFDP(2,4)=0
       DFDP(3,4)=0

       DFDP(1,5)=x**3/r
       DFDP(2,5)=0
       DFDP(3,5)=0

       DFDP(1,6)=(y-x)**3 / r
       DFDP(2,6)=-(y-x)**3
       DFDP(3,6)=0

      END SUBROUTINE FUNC

      SUBROUTINE STPNT(NDIM,U,PAR,T)
!     ---------- ----- 

      IMPLICIT NONE
      INTEGER, INTENT(IN) :: NDIM
      DOUBLE PRECISION, INTENT(INOUT) :: U(NDIM),PAR(*)
      DOUBLE PRECISION, INTENT(IN) :: T

!----------------------------------------------------------------------
! Problem parameters (only PAR(1-9) are available to the user) :

       PAR(1)=-0.721309D0         ! nu
       PAR(2)=0.6                 ! beta
       PAR(3)=0.0                 ! gamma
       PAR(4)=0.6                 ! r
       PAR(5)=0.328578            ! a_3
       PAR(6)=0.933578            ! b_3

       PAR(11)=36.13              ! Truncated time interval 

!----------------------------------------------------------------------
! If IEQUIB >0 put initial equilibrium in PAR(11+i), i=1,...,NDIM :

      PAR(12) = 0.0
      PAR(13) = 0.0
      PAR(14) = 0.0

      END SUBROUTINE STPNT

      SUBROUTINE PVLS(NDIM,U,PAR)
!     ---------- ----

      IMPLICIT NONE
      INTEGER, INTENT(IN) :: NDIM
      DOUBLE PRECISION, INTENT(IN) :: U(NDIM)
      DOUBLE PRECISION, INTENT(INOUT) :: PAR(*)
! Homoclinic bifurcations COMMON block needed here :
      INTEGER ITWIST,ISTART,IEQUIB,NFIXED,NPSI,NUNSTAB,NSTAB,NREV
      COMMON /BLHOM/ ITWIST,ISTART,IEQUIB,NFIXED,NPSI,NUNSTAB,NSTAB,NREV

      INTEGER I

! If IEQUIB =0 put analytic equilibrium in PAR(11+i), i=1..NDIM

      IF(IEQUIB.EQ.0)THEN
        DO I=1,NDIM
          PAR(11+I)= 0.0
        ENDDO
      ENDIF

      END SUBROUTINE PVLS

      SUBROUTINE BCND 
      END SUBROUTINE BCND

      SUBROUTINE ICND 
      END SUBROUTINE ICND

      SUBROUTINE FOPT 
      END SUBROUTINE FOPT