[updated 10.Jul.2005]
Librairie ef3d > Fichier qulkfb.f |
SUBROUTINE QULKFB (K,KELFCT,X1,X2,P,DP,DDP)
Auteur : D.Martin & G. Vial (Septembre 2000)
Refonte complete pour les elements de Lagrange aux abcisses quelconques
(en particulier de Gauss-Lobatto) : D.Martin (29 Aout 2001)
Derniere modification : D.Martin (6 Septembre 2001)
Version 1.1.1
Fonctions de base et derivees - Quadrangle de Lagrange d'ordre k
construites par produit tensoriel a partir des fonctions de base 1D
-- Arguments d'entree --
K ordre de l'element de Lagrange
KELFCT variable CHARACTER indiquant
si KELFCT(1:1).NE.'N' Calcul des fonctions de base
si KELFCT(2:2).NE.'N' Calcul des derivees 1eres des fonctions de base
si KELFCT(3:3).NE.'N' Calcul des derivees 2ndes des fonctions de base
X1,X2 coordonnees du point ou on demande le calcul
-- Arguments de sortie --
P tableau contenant les valeurs des fonctions au point X1,X2
DP - - - - derivees premieres
DDP - - - - derivees secondes
-- Numerotations locales --
Segments
2----1 2---3---1 2---4---3---1 2---5---4---3---1 2---6---5---4---3---1
k=1 k=2 k=3 k=4 k=5
Quadrangles
3----2 3---6---2 3--10---6---2 3--14--10---6---2 3--18--14--10---6---2
| | | | | | | | | |
4----1 7 9 5 7 15 14 9 7 19 22 18 13 7 23 30 26 22 17
k=1 | | | | | | | |
4---8---1 11 16 13 5 11 23 25 21 9 11 27 35 34 29 13
k=2 | | | | | |
4---8--12---1 15 20 24 17 5 15 31 36 33 25 9
k=3 | | | |
4---8--12--16---1 19 24 28 32 21 5
k=4 | |
4---8--12--16--20---1
k=5
IMPLICIT NONE
CHARACTER*(*) KELFCT
INTEGER K
REAL X1,X2,P(*),DP(2,*),DDP(3,*)
INTEGER LAGDMX
PARAMETER ( LAGDMX = 100)
REAL P1(LAGDMX+1),DP1(LAGDMX+1),DDP1
REAL P2(LAGDMX+1),DP2(LAGDMX+1),DDP2
INTEGER KNEW,Q,NK,NE,I1,I2,I3,I33,IK,IKK,I,J
Fonctions de l'element 1D
IF (KELFCT(1:2).EQ.'YN') THEN
CALL SELKFB (K,'YNN',X1,P1,DP1,DDP1)
CALL SELKFB (K,'YNN',X2,P2,DP2,DDP2)
ELSE
CALL SELKFB (K,'YYN',X1,P1,DP1,DDP1)
CALL SELKFB (K,'YYN',X2,P2,DP2,DDP2)
ENDIF
KNEW=K
IF (KNEW.GT.LAGDMX) KNEW=KNEW-LAGDMX
Q=(KNEW/2)+1
Fonctions de base
IF (KELFCT(1:1).NE.'N') THEN
NK=KNEW
NE=0
I1=1
I2=2
I3=3
IK=NK+1
DO 12 I=1,Q
IF (NK.GT.0) THEN
Les 4 sommets du 'perimetre' en cours
NE=NE+1
P(NE)=P1(I1)*P2(I2)
NE=NE+1
P(NE)=P1(I1)*P2(I1)
NE=NE+1
P(NE)=P1(I2)*P2(I1)
NE=NE+1
P(NE)=P1(I2)*P2(I2)
Les noeuds 'internes' aux aretes du 'perimetre' en cours
I33=I3
IKK=IK
DO 11 J=2,NK
NE=NE+1
P(NE)=P1(I1)*P2(IKK)
NE=NE+1
P(NE)=P1(I33)*P2(I1)
NE=NE+1
P(NE)=P1(I2)*P2(I33)
NE=NE+1
P(NE)=P1(IKK)*P2(I2)
IKK=IKK-1
I33=I33+1
11 CONTINUE
I3=I3+1
IK=IK-1
I1=I1+1
I2=I2-1
IF (I.EQ.1) THEN
I1=3
I2=NK+1
ENDIF
NK=NK-2
ELSE
Le point central le cas echeant
NE=NE+1
P(NE)=P1(I1)*P2(I1)
ENDIF
12 CONTINUE
ENDIF
Derivees premieres
IF (KELFCT(2:2).NE.'N') THEN
NK=KNEW
NE=0
I1=1
I2=2
I3=3
IK=NK+1
DO 22 I=1,Q
IF (NK.GT.0) THEN
Les 4 sommets du 'perimetre' en cours
NE=NE+1
DP(1,NE)=DP1(I1)*P2(I2)
DP(2,NE)=P1(I1)*DP2(I2)
NE=NE+1
DP(1,NE)=DP1(I1)*P2(I1)
DP(2,NE)=P1(I1)*DP2(I1)
NE=NE+1
DP(1,NE)=DP1(I2)*P2(I1)
DP(2,NE)=P1(I2)*DP2(I1)
NE=NE+1
DP(1,NE)=DP1(I2)*P2(I2)
DP(2,NE)=P1(I2)*DP2(I2)
I33=I3
IKK=IK
DO 21 J=2,NK
NE=NE+1
DP(1,NE)=DP1(I1)*P2(IKK)
DP(2,NE)=P1(I1)*DP2(IKK)
NE=NE+1
DP(1,NE)=DP1(I33)*P2(I1)
DP(2,NE)=P1(I33)*DP2(I1)
NE=NE+1
DP(1,NE)=DP1(I2)*P2(I33)
DP(2,NE)=P1(I2)*DP2(I33)
NE=NE+1
DP(1,NE)=DP1(IKK)*P2(I2)
DP(2,NE)=P1(IKK)*DP2(I2)
IKK=IKK-1
I33=I33+1
21 CONTINUE
I3=I3+1
IK=IK-1
I1=I1+1
I2=I2-1
IF (I.EQ.1) THEN
I1=3
I2=NK+1
ENDIF
NK=NK-2
ELSE
Le point central le cas echeant
NE=NE+1
DP(1,NE)=DP1(I1)*P2(I1)
DP(2,NE)=P1(I1)*DP2(I1)
ENDIF
22 CONTINUE
ENDIF
END
qulkfb est appelé dans