1*59599516SKenneth E. Jansen /*-------------------------------------------------------------------------
2*59599516SKenneth E. Jansen Scientific Computation Research Center, RPI, Troy NY
3*59599516SKenneth E. Jansen (C) Copyright 1995, RPI-SCOREC
4*59599516SKenneth E. Jansen
5*59599516SKenneth E. Jansen Project : shapeFuntions
6*59599516SKenneth E. Jansen Author(s): Saikat Dey
7*59599516SKenneth E. Jansen Creation : Oct., 95
8*59599516SKenneth E. Jansen Modifi. :
9*59599516SKenneth E. Jansen Function :
10*59599516SKenneth E. Jansen return derivative hierarchic mode shapes associated with
11*59599516SKenneth E. Jansen a given mesh entity of a specified polynomial order.
12*59599516SKenneth E. Jansen -------------------------------------------------------------------------*/
13*59599516SKenneth E. Jansen
14*59599516SKenneth E. Jansen #include <math.h>
15*59599516SKenneth E. Jansen #include "shapeFuncInternals.h"
16*59599516SKenneth E. Jansen
17*59599516SKenneth E. Jansen #ifdef __cplusplus
18*59599516SKenneth E. Jansen extern "C" {
19*59599516SKenneth E. Jansen #endif
20*59599516SKenneth E. Jansen
E_modeShapeDrv(int p,double * L,double drv[2])21*59599516SKenneth E. Jansen int E_modeShapeDrv(int p, double *L, double drv[2]) {
22*59599516SKenneth E. Jansen /* Return derivative edge mode shape function evaluated along an edge at
23*59599516SKenneth E. Jansen a point for spectral interpolation order p L[0,1] in [0,1] r+s <= 1
24*59599516SKenneth E. Jansen */
25*59599516SKenneth E. Jansen return EnDrv(p-2,L[0],L[1],drv);
26*59599516SKenneth E. Jansen }
27*59599516SKenneth E. Jansen
F_modeShapeTriDrv(int p,int i,double * L,double mdrv[2])28*59599516SKenneth E. Jansen int F_modeShapeTriDrv(int p, int i, double *L, double mdrv[2]) {
29*59599516SKenneth E. Jansen int alpha,beta,found,count;
30*59599516SKenneth E. Jansen double rs,rst2,P1P2,t,P1,P2,dP1,dP2,t1,t2;
31*59599516SKenneth E. Jansen /* return i-th triangular face mode derivative of polynomial order p
32*59599516SKenneth E. Jansen note: there are p-2 modes with polynomial order p */
33*59599516SKenneth E. Jansen
34*59599516SKenneth E. Jansen if( p < 3 || i < 0 || i > p-3 )
35*59599516SKenneth E. Jansen return 0.0 ;
36*59599516SKenneth E. Jansen
37*59599516SKenneth E. Jansen count = found = 0;
38*59599516SKenneth E. Jansen for(alpha=0; alpha <= p-3; alpha++) {
39*59599516SKenneth E. Jansen for(beta=0; beta <= p-3; beta++) {
40*59599516SKenneth E. Jansen if( alpha+beta == p-3 ) {
41*59599516SKenneth E. Jansen if( count == i )
42*59599516SKenneth E. Jansen found=1;
43*59599516SKenneth E. Jansen else
44*59599516SKenneth E. Jansen count++;
45*59599516SKenneth E. Jansen }
46*59599516SKenneth E. Jansen if(found)
47*59599516SKenneth E. Jansen break ;
48*59599516SKenneth E. Jansen }
49*59599516SKenneth E. Jansen if(found)
50*59599516SKenneth E. Jansen break;
51*59599516SKenneth E. Jansen }
52*59599516SKenneth E. Jansen if( found ) {
53*59599516SKenneth E. Jansen return FnDrv(alpha,beta,L[0],L[1],mdrv);
54*59599516SKenneth E. Jansen } else
55*59599516SKenneth E. Jansen return 0;
56*59599516SKenneth E. Jansen }
57*59599516SKenneth E. Jansen
F_modeShapeQuadDrv(int p,int i,double * L,double mdrv[2])58*59599516SKenneth E. Jansen int F_modeShapeQuadDrv(int p, int i, double *L, double mdrv[2]) {
59*59599516SKenneth E. Jansen return 0;
60*59599516SKenneth E. Jansen }
61*59599516SKenneth E. Jansen
R_modeShapeTetDrv(int p,int i,double * L,double mdrv[3])62*59599516SKenneth E. Jansen int R_modeShapeTetDrv(int p, int i, double *L, double mdrv[3]) {
63*59599516SKenneth E. Jansen int alpha,beta,gamma,count,found ;
64*59599516SKenneth E. Jansen double Pr,Ps,Pt,dPr,dPs,dPt,w,rst,rstw2,PrPsPt,t1;
65*59599516SKenneth E. Jansen
66*59599516SKenneth E. Jansen /* return the i-th mode shape of polynomial order p for this region , there are
67*59599516SKenneth E. Jansen (p-2)*(p-3)/2 mode shapes of polynomial order p */
68*59599516SKenneth E. Jansen if( p < 4 || i < 0 || i > (((p-2)*(p-3)/2)-1) )
69*59599516SKenneth E. Jansen return 0.0;
70*59599516SKenneth E. Jansen
71*59599516SKenneth E. Jansen count = 0;
72*59599516SKenneth E. Jansen found = 0;
73*59599516SKenneth E. Jansen for( alpha=0; alpha <= p-4; alpha++) {
74*59599516SKenneth E. Jansen for( beta=0; beta <= p-4; beta++) {
75*59599516SKenneth E. Jansen for( gamma=0; gamma <= p-4; gamma++) {
76*59599516SKenneth E. Jansen if( alpha+beta+gamma == p-4 ) {
77*59599516SKenneth E. Jansen if( count == i )
78*59599516SKenneth E. Jansen found = 1;
79*59599516SKenneth E. Jansen else
80*59599516SKenneth E. Jansen count++;
81*59599516SKenneth E. Jansen }
82*59599516SKenneth E. Jansen if(found) break;
83*59599516SKenneth E. Jansen }
84*59599516SKenneth E. Jansen if(found) break;
85*59599516SKenneth E. Jansen }
86*59599516SKenneth E. Jansen if(found) break;
87*59599516SKenneth E. Jansen }
88*59599516SKenneth E. Jansen if( found ) {
89*59599516SKenneth E. Jansen return BnDrv(alpha,beta,gamma,L[0],L[1],L[2],mdrv);
90*59599516SKenneth E. Jansen } else
91*59599516SKenneth E. Jansen return 0;
92*59599516SKenneth E. Jansen }
93*59599516SKenneth E. Jansen
R_modeShapeHexDrv(int p,int i,double * L,double mdrv[3])94*59599516SKenneth E. Jansen int R_modeShapeHexDrv(int p, int i, double *L, double mdrv[3]) {
95*59599516SKenneth E. Jansen return 0;
96*59599516SKenneth E. Jansen }
97*59599516SKenneth E. Jansen
98*59599516SKenneth E. Jansen #ifdef __cplusplus
99*59599516SKenneth E. Jansen }
100*59599516SKenneth E. Jansen #endif
101*59599516SKenneth E. Jansen
102*59599516SKenneth E. Jansen
103*59599516SKenneth E. Jansen
104*59599516SKenneth E. Jansen
105*59599516SKenneth E. Jansen
106*59599516SKenneth E. Jansen
107