c Rutina 053
c-------------------------------------------------------------CALCONS---
C
SUBROUTINE CALCONS(NTIP,E,POIS,G12,B11,B22,B12,B66,GAMMA,ACU,BMAS,
@ K11,K22,K12,C11,C22,A1,A2,B1,B2,C1,C2,D1,D2,
@ BMEN,AMAS,AMEN,MU1,MU2,P,Q,CUK,QP,PP,
@ QQP,QQQ,UU1,UU2,TT1,TT2,DD1,DD2,SS1,SS2)
C
C Esta rutina calcula las constantes que aparecen en la ley de com_
C portamiento y en la solucion fundamental de materiales ortotropos.
C
IMPLICIT NONE
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C VARIABLES EXTERNAS (ENTRADA)
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C NTIP ................... tipo de problema:
C 0 ==> deformacion plana
C 1 ==> tension plana
INTEGER NTIP
C E(3) ........... modulos de elasticidad E11, E22, E33
C POIS(3) ........ coeficiente de Poison POIS12, POIS13, POIS23
C G12 ........... modulo de cortadura G12
REAL*8 E(3),POIS(3),G12
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C VARIABLES EXTERNAS (SALIDA)
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C GAMMA() ................ constante c (solucion fundamental)
REAL*8 GAMMA
C BMAS() ................. constante b (solucion fundamental)
REAL*8 BMAS
C B11 ........... coeficiente B11 de la ley de comportamiento
C B12 ........... coeficiente B12 de la ley de comportamiento
C B22 ........... coeficiente B22 de la ley de comportamiento
C B66 ........... coeficiente B66 de la ley de comportamiento
C
REAL*8 B11,B12,B22,B66
C ACU ........... constante a^2 (solucion fundamental)
REAL*8 ACU
C K11 ........... constantes K11 (solucion fundamental)
C K22 ........... constantes K22 (solucion fundamental)
C K12 ........... constantes K12 (solucion fundamental)
C C11 ........... constantes C11 (solucion fundamental)
C C22 ........... constantes C22 (solucion fundamental)
REAL*8 K11,K22,K12
REAL*8 C11,C22
C A1 ............ constantes A1 (solucion fundamental)
C A2 ............ constantes A2 (solucion fundamental)
C B1 ............ constantes B1 (solucion fundamental)
C B2 ............ constantes B2 (solucion fundamental)
C C1 ............ constantes C1 (solucion fundamental)
C C2 ............ constantes C2 (solucion fundamental)
C D1 ............ constantes D1 (solucion fundamental)
C D2 ............ constantes D2 (solucion fundamental)
REAL*8 A1,A2,B1,B2
REAL*8 C1,C2,D1,D2
REAL*8 AMAS,AMEN
complex*16 P(2,2),Q(2,2),CUK(2)
complex*16 BMEN,MU1,MU2
complex*16 QP(2,2,2),PP(2,2,2)
complex*16 QQQ(3,2,2),QQP(3,2,2)
complex*16 UU1(2,2),UU2(2,2),TT1(2,2)
complex*16 TT2(2,2),DD1(3,2),DD2(3,2)
complex*16 SS1(3,2),SS2(3,2)
C
CCCCCCCCCCCCCCCCCCCCCCCCCCC
C VARIABLES INTERNAS
CCCCCCCCCCCCCCCCCCCCCCCCCCC
C
INTEGER I,J,K,L,N
REAL*8 PI,AUX1,AUX2,AUX3,RACU
COMPLEX*16 IC
C PI .... sin comentario
C AUX1,AUX2,AUX3 vbles auxiliares para definir las constantes
C IC .... unidad imaginaria
IC=(0.D00,1.D00)
PI=4.D00*DATAN(1.D00)
C
CCCCCCCCCCCCCCC
C RUTINA
CCCCCCCCCCCCCCC
C
IF (NTIP.EQ.1) THEN ! tension plana
B11=1/E(1)
B22=1/E(2)
B12=-POIS(1)/E(1)
B66=1/G12
ELSE ! deformacion plana
B11=1/E(1)-POIS(2)**2*E(3)/E(1)**2
B22=1/E(2)-POIS(3)**2*E(3)/E(2)**2
B12=-POIS(1)/E(1)-POIS(2)*POIS(3)*E(3)/
* E(1)/E(2)
B66=1/G12
ENDIF
GAMMA=B66/2
ACU=DSQRT(B11*B22)-B12-GAMMA
BMAS=DSQRT(DSQRT(B11*B22)+B12+GAMMA)
IF (DABS(ACU).LT.1D-10) ACU=0
AUX1=8.D00*DSQRT(2.D00)*PI*BMAS*DSQRT(B22)
AUX2=8.D00*DSQRT(2.D00)*PI*BMAS*DSQRT(B11)
AUX3=(BMAS**2-GAMMA)*GAMMA+(ACU+GAMMA)*BMAS**2
K11=AUX3/AUX1
K22=AUX3/AUX2
IF (ACU.GT.0.D00) THEN ! a real
AUX1=4.D00*DSQRT(2.D00)*PI*DSQRT(ACU*B22)
AUX2=4.D00*DSQRT(2.D00)*PI*DSQRT(ACU*B11)
AUX3=(BMAS**2-GAMMA)*ACU-(ACU+GAMMA)*GAMMA
C11=AUX3/AUX1
C22=AUX3/AUX2
K12=(ACU*BMAS**2+GAMMA**2)/(8.D00*PI*DSQRT(ACU)*BMAS)
ELSE IF (ACU.LT.0.D00) THEN !a imaginario puro
AUX1=8.D00*DSQRT(2.D00)*PI*DSQRT(-ACU*B22)
AUX2=8.D00*DSQRT(2.D00)*PI*DSQRT(-ACU*B11)
AUX3=(BMAS**2-GAMMA)*ACU-(ACU+GAMMA)*GAMMA
C11=AUX3/AUX1
C22=AUX3/AUX2
K12=(ACU*BMAS**2+GAMMA**2)/(4.D0*PI*DSQRT(-ACU)*BMAS)
ELSE
AUX1=4.D00*DSQRT(2.D00)*PI*DSQRT(B22)
AUX2=4.D00*DSQRT(2.D00)*PI*DSQRT(B11)
C11=GAMMA**2/AUX1
C22=GAMMA**2/AUX2
K12=GAMMA**2/(4.D00*PI*BMAS)
ENDIF
A1=DSQRT(B22)*(BMAS**2-GAMMA)/(2.D00*DSQRT(2.D00)*PI*BMAS)
A2=DSQRT(B11)*(BMAS**2-GAMMA)/(2.D00*DSQRT(2.D00)*PI*BMAS)
B1=DSQRT(B11)*(BMAS**2+GAMMA)/(2.D00*DSQRT(2.D00)*PI*BMAS)
B2=DSQRT(B22)*(BMAS**2+GAMMA)/(2.D00*DSQRT(2.D00)*PI*BMAS)
C1=((BMAS**2+GAMMA)*DSQRT(B11*B22)-2.D00*(ACU+GAMMA)*
* BMAS**2)/(2.D00*DSQRT(2.D00)*PI*BMAS*DSQRT(B11))
C2=(BMAS**2-GAMMA)*B11/(2.D00*DSQRT(2.D00)*PI*BMAS*
* DSQRT(B22))
D1=(BMAS**2-GAMMA)*B22/(2.D00*DSQRT(2.D00)*PI*BMAS*
* DSQRT(B11))
D2=((BMAS**2+GAMMA)*DSQRT(B11*B22)-2.D00*(ACU+GAMMA)*
* BMAS**2)/(2.D00*DSQRT(2.D00)*PI*BMAS*DSQRT(B22))
IF(ACU.NE.0.D00) THEN
IF(ACU.LT.0.D00) THEN
RACU=DSQRT(-ACU)
BMEN=IC*RACU
ELSE
RACU=DSQRT(ACU)
BMEN=RACU*(1.D00,0.D00)
END IF
ELSE
BMEN=(0.D00,0.D00)
END IF
AMAS=DSQRT(B11*B22)+B12
AMEN=DSQRT(B11*B22)-B12
IF(ACU.NE.0.D00) THEN
MU1=(BMEN+IC*BMAS)/DSQRT(2.D00*B11)
MU2=(-BMEN+IC*BMAS)/DSQRT(2.D00*B11)
ELSE
MU1=IC*DSQRT(DSQRT(B22/B11))
MU2=MU1
END IF
IF(.NOT.(MU1.EQ.MU2)) THEN ! raices del pol. caract. <>
P(1,1)=-GAMMA+IC*BMEN*BMAS
P(1,2)=-GAMMA-IC*BMEN*BMAS
P(2,1)=-MU1*(GAMMA+IC*BMEN*BMAS)
P(2,2)=-MU2*(GAMMA-IC*BMEN*BMAS)
Q(1,1)=-MU1
Q(1,2)=-MU2
Q(2,1)=(1.D00,0.D00)
Q(2,2)=(1.D00,0.D00)
CUK(1)=-1.D00/(4.D00*IC*BMEN*BMAS*MU1)
CUK(2)=1.D00/(4.D00*IC*BMEN*BMAS*MU2)
DO I=1,2
DO J=1,2
DO L=1,2
QP(J,I,L)=Q(J,L)*P(I,L)*CUK(L)/(PI*IC)
PP(I,J,L)=P(I,L)*P(J,L)*CUK(L)/(PI*IC)
END DO
END DO
END DO
DO N=1,2 ! N es la coordenada segun los ejes 1 y 2
DO L=1,2 ! L representa la raiz del pol. caract.
QQP(1,N,L)=Q(1,L)*Q(1,L)*P(N,L)*CUK(L)/(PI*IC)
QQP(2,N,L)=Q(1,L)*Q(2,L)*P(N,L)*CUK(L)/(PI*IC)
QQP(3,N,L)=Q(2,L)*Q(2,L)*P(N,L)*CUK(L)/(PI*IC)
QQQ(1,N,L)=Q(1,L)*Q(1,L)*Q(N,L)*CUK(L)/(PI*IC)
QQQ(2,N,L)=Q(1,L)*Q(2,L)*Q(N,L)*CUK(L)/(PI*IC)
QQQ(3,N,L)=Q(2,L)*Q(2,L)*Q(N,L)*CUK(L)/(PI*IC)
END DO
END DO
ELSE ! si BETA MENOS = 0 se calculan los nuevos nucleos
DO K=1,2
DO L=1,2
IF(K.EQ.L) THEN
UU1(K,L)=((-1.D00)**(K))*(MU1**(K+L-3))*GAMMA*
* (GAMMA-2*BMAS*BMAS)/(2.D00*PI*IC*
* BMAS*BMAS)
ELSE
UU1(K,L)=(0.D00,0.D00)
END IF
UU2(K,L)=-(MU1**(K+L-3))*GAMMA*GAMMA/(4.D00*
* PI*IC*BMAS*BMAS)
TT1(K,L)=IC*(((-AMAS)/(BMAS*BMAS))**
* ABS(K-L))*(MU1**(K-L))/(2.D00*PI)
TT2(K,L)=((-1.D00)**(L))*IC*GAMMA*(MU1**(K-L))/
* (4.D00*PI*BMAS*BMAS)
END DO
END DO
C Nucleos para el calculo de las tensiones en puntos internos
DO J=1,3
DO L=1,2
DD1(J,L)=((-1.D00)**(J+1))*(((0.5D00-J+L)*GAMMA/
* (BMAS*BMAS))-((-1.D00)**(L)))*
* ((MU1)**(1-J+L))/(2.D00*PI*IC)
DD2(J,L)=((-1.D00)**(J+1))*GAMMA*((MU1)**
* (1-J+L))/(4.D00*PI*IC*BMAS*BMAS)
SS1(J,L)=((-1.D00)**(J+L+1))*(3.5D00-J-L)*((MU1)**
* (4-J-L))/(2.D00*PI*IC*BMAS*BMAS)
SS2(J,L)=((-1.D00)**(J+L+1))*((MU1)**(4-J-L))/
* (4.D00*PI*IC*BMAS*BMAS)
END DO
END DO
END IF
C
RETURN
C
END
C
e-REdING. Biblioteca de la Escuela Superior de Ingenieros de Sevilla.
