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 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189
|
SUBROUTINE DPSE4
C
C PRESSURE STIFFNESS CALCULATIONS FOR A QUADRILATERAL MEMBRANE
C ELEMENT, WHICH HAS 4 GRID POINTS.
C THREE 6X6 STIFFNESS MATRICES FOR THE PIVOT POINT ARE INSERTED.
C
C DOUBLE PRECISION VERSION
C
C WRITTEN BY E. R. CHRISTENSEN/SVERDRUP, 9/91, VERSION 1.1
C INSTALLED IN NASTRAN AS ELEMENT DPSE4 BY G.CHAN/UNISYS, 2/92
C
C REFERENCE - E. CHRISTENEN: 'ADVACED SOLID ROCKET MOTOR (ASRM)
C MATH MODELS - PRESSURE STIFFNESS EFFECTS ANALYSIS',
C NASA TD 612-001-02, AUGUST 1991
C
C LIMITATION -
C (1) ALL GRID POINTS USED BY ANY OF THE CPSE2/3/4 ELEMENTS MUST BE
C IN BASIC COORDINATE SYSTEM!!!
C (2) CONSTANT PRESSURE APPLIED OVER AN ENCLOSED VOLUMN ENCOMPASSED
C BY THE CPSE2/3/4 ELEMENTRS
C (3) PRESSURE ACTS NORMALLY TO THE CPSE2/3/4 SURFACES
C
C SEE NASTRAN DEMONSTRATION PROBLEM - T13022A
C
DOUBLE PRECISION GAMMA,KIJ,DP,C,SIGN
DIMENSION NECPT(6)
C COMMON /SYSTEM/ IBUF,NOUT
COMMON /DS1AAA/ NPVT,ICSTM,NCSTM
COMMON /DS1AET/ ECPT(26),DUM2(2),DUM12(12)
COMMON /DS1ADP/ GAMMA,KIJ(36),DP(26),C(12),SIGN(3),NK(3),IK(3)
EQUIVALENCE (NECPT(1),ECPT(1))
C
C ECPT FOR THE PRESSURE STIFFNESS CPES4 ELEMENT
C
C ECPT( 1) = ELEMENT ID
C ECPT( 2) = SIL FOR GRID POINT A OR 1
C ECPT( 3) = SIL FOR GRID POINT B OR 2
C ECPT( 4) = SIL FOR GRID POINT C OR 3
C ECPT( 5) = SIL FOR GRID POINT C OR 4
C ECPT( 6) = PRESSURE
C ECPT( 7) = NOT USED
C ECPT( 8) = NOT USED
C ECPT( 9) = NOT USED
C ECPT(10) = COORD. SYSTEM ID 1
C ECPT(11) = X1
C ECPT(12) = Y1
C ECPT(13) = Z1
C ECPT(14) = COORD. SYSTEM ID 2
C ECPT(15) = X2
C ECPT(16) = Y2
C ECPT(17) = Z2
C ECPT(18) = COORD. SYSTEM ID 3
C ECPT(19) = X3
C ECPT(20) = Y3
C ECPT(21) = Z3
C ECPT(22) = COORD. SYSTEM ID 4
C ECPT(23) = X4
C ECPT(24) = Y4
C ECPT(25) = Z4
C ECPT(26) = ELEMENT TEMPERATURE
C ECPT(27) THRU ECPT(40) = DUM2 AND DUM12, NOT USED IN THIS ROUTINE
C
C STORE ECPT IN DOUBLE PRECISION
C
DP(6) = ECPT(6)
K = 10
DO 20 I = 1,4
DO 10 J = 1,3
K = K + 1
10 DP(K) = ECPT(K)
20 K = K + 1
C
C CALCULATE THE FOUR VECTORS GAB, GAC, GAD, AND GBD USED IN
C COMPUTING THE PRESSURE STIFFNESS MATRIC
C
C GAB = RA + RB - RC - RD
C GAC = RB - RD
C GAD =-RA + RB + RC - RD
C GBD =-RA + RC
C
C GAB STORED IN C( 1), C( 2), C( 3)
C GAC STORED IN C( 4), C( 5), C( 6)
C GAD STORED IN C( 7), C( 8), C( 9)
C GBD STORED IN C(10), C(11), C(12)
C
C(1) = DP(11) + DP(15) - DP(19) - DP(23)
C(2) = DP(12) + DP(16) - DP(20) - DP(24)
C(3) = DP(13) + DP(17) - DP(21) - DP(25)
C
C(4) = DP(15) - DP(23)
C(5) = DP(16) - DP(24)
C(6) = DP(17) - DP(25)
C
C(7) =-DP(11) + DP(15) + DP(19) - DP(23)
C(8) =-DP(12) + DP(16) + DP(20) - DP(24)
C(9) =-DP(13) + DP(17) + DP(21) - DP(25)
C
C(10)=-DP(11) + DP(19)
C(11)=-DP(12) + DP(20)
C(12)=-DP(13) + DP(21)
C
DO 30 I = 1,4
IF (NECPT(I+1) .NE. NPVT) GO TO 30
NPIVOT = I
GO TO 40
30 CONTINUE
RETURN
C
C GENERATE THE THREE BY THREE PARTITIONS IN GLOBAL COORDINATES HERE
C
C SET COUNTERS ACCORDING TO WHICH GRID POINT IS THE PIVOT
C
40 IF (NPIVOT .EQ. 4) GO TO 80
IF (NPIVOT-2) 50,60,70
C
C SET COUNTERS AND POINTERS FOR CALCULATING KAB, KAC, KAD
C
50 NK(1) = 2
NK(2) = 3
NK(3) = 4
IK(1) = 1
IK(2) = 4
IK(3) = 7
SIGN(1) = 1.0D0
SIGN(2) = 1.0D0
SIGN(3) = 1.0D0
GO TO 90
C
C SET COUNTERS AND POINTERS FOR CALCULATING KBA, KBC, KBD
C NOTE THAT KBA = -KAB
C
60 NK(1) = 1
NK(2) = 3
NK(3) = 4
IK(1) = 1
IK(2) = 7
IK(3) = 10
SIGN(1) =-1.0D0
SIGN(2) = 1.0D0
SIGN(3) = 1.0D0
GO TO 90
C
C SET COUNTERS AND POINTERS FOR CALCULATING KCA, KCB, KCD
C NOTE THAT KCA = -KAC, KCB = -KBC
C
70 NK(1) = 1
NK(2) = 2
NK(3) = 4
IK(1) = 4
IK(2) = 7
IK(3) = 1
SIGN(1) =-1.0D0
SIGN(2) =-1.0D0
SIGN(3) =-1.0D0
GO TO 90
C
80 NK(1) = 1
NK(2) = 2
NK(3) = 3
IK(1) = 7
IK(2) = 10
IK(3) = 1
SIGN(1) =-1.0D0
SIGN(2) =-1.0D0
SIGN(3) = 1.0D0
C
90 GAMMA =-DP(6)/12.0D0
DO 110 I = 1,3
DO 100 J = 1,36
100 KIJ(J) = 0.0D0
K1 = IK(I)
K2 = K1 + 1
K3 = K1 + 2
SG = GAMMA*SIGN(I)
KIJ( 2) =-C(K3)*SG
KIJ( 3) = C(K2)*SG
KIJ( 7) = C(K3)*SG
KIJ( 9) =-C(K1)*SG
KIJ(13) =-C(K2)*SG
KIJ(14) = C(K1)*SG
C
C ASSEMBLE INTO THE GLOBAL STIFFNESS MATRIX
C
IAS = NK(I)
CALL DS1B (KIJ(1),NECPT(IAS+1))
110 CONTINUE
C
RETURN
END
|