320 INTEGER(INTG),
INTENT(OUT) :: ERR
324 enters(
"BASES_FINALISE",err,error,*999)
334 exits(
"BASES_FINALISE")
336 999 errorsexits(
"BASES_FINALISE",err,error)
349 INTEGER(INTG),
INTENT(OUT) :: ERR
353 enters(
"BASES_INITIALISE",err,error,*999)
358 exits(
"BASES_INITIALISE")
360 999 errorsexits(
"BASES_INITIALISE",err,error)
374 INTEGER(INTG),
INTENT(OUT) :: ERR
377 INTEGER(INTG) :: ni,nic,nn,nn1,nn2,nn3,nn4,ns,local_line_idx,local_face_idx
380 enters(
"BASIS_CREATE_FINISH",err,error,*999)
382 IF(
ASSOCIATED(basis))
THEN 383 SELECT CASE(basis%TYPE)
391 local_error=
"Basis type "//
trim(
number_to_vstring(basis%TYPE,
"*",err,error))//
" is invalid or not implemented" 392 CALL flagerror(local_error,err,error,*999)
394 basis%BASIS_FINISHED=.true.
396 CALL flagerror(
"Basis is not associated",err,error,*999)
409 &
'(" Interpolation type(nic):",4(X,I2))',
'(25X,4(X,I2))',err,error,*999)
411 &
'(" Interpolation order(nic):",4(X,I2))',
'(26X,4(X,I2))',err,error,*999)
413 &
'(" Collapsed xi(ni):",3(X,I2))',
'(26X,3(X,I2))',err,error,*999)
418 &
'(" Number of nodes(nic):",4(X,I2))',
'(22X,4(X,I2))',err,error,*999)
420 &
'(" Number of derivatives(nn):",8(X,I2))',
'(28X,8(X,I2))',err,error,*999)
425 &
'(" Node at collapse(nn):",8(X,L1))',
'(23X,8(X,L1))',err,error,*999)
427 DO nic=1,basis%NUMBER_OF_XI_COORDINATES
430 &
'(" INDEX(nn) :",16(X,I2))',
'(18X,16(X,I2))',err,error,*999)
434 SELECT CASE(basis%NUMBER_OF_XI_COORDINATES)
436 DO nn1=1,basis%NUMBER_OF_NODES_XIC(1)
442 DO nn2=1,basis%NUMBER_OF_NODES_XIC(2)
444 DO nn1=1,basis%NUMBER_OF_NODES_XIC(1)
451 DO nn3=1,basis%NUMBER_OF_NODES_XIC(3)
453 DO nn2=1,basis%NUMBER_OF_NODES_XIC(2)
455 DO nn1=1,basis%NUMBER_OF_NODES_XIC(1)
463 DO nn4=1,basis%NUMBER_OF_NODES_XIC(4)
465 DO nn3=1,basis%NUMBER_OF_NODES_XIC(3)
467 DO nn2=1,basis%NUMBER_OF_NODES_XIC(2)
469 DO nn1=1,basis%NUMBER_OF_NODES_XIC(1)
472 & ,basis%NODE_POSITION_INDEX_INV(nn1,nn2,nn3,nn4),err,error,*999)
478 CALL flagerror(
"Invalid number of xi coordinates",err,error,*999)
481 DO ni=1,basis%NUMBER_OF_XI
483 DO nn=1,basis%NUMBER_OF_NODES
486 & basis%DERIVATIVE_ORDER_INDEX(:,nn,ni),
'(" INDEX(nk):",8(X,I2))',
'(18X,8(X,I2))',err,error,*999)
491 DO nn=1,basis%NUMBER_OF_NODES
494 & basis%ELEMENT_PARAMETER_INDEX(:,nn),
'(" INDEX(nk) :",8(X,I2))',
'(18X,8(X,I2))',err,error,*999)
497 DO ns=1,basis%NUMBER_OF_ELEMENT_PARAMETERS
500 & basis%ELEMENT_PARAMETER_INDEX_INV(:,ns),
'(" INDEX(:) :",2(X,I2))',
'(18X,2(X,I2))',err,error,*999)
503 DO nn=1,basis%NUMBER_OF_NODES
506 & basis%PARTIAL_DERIVATIVE_INDEX(:,nn),
'(" INDEX(nk) :",8(X,I2))',
'(18X,8(X,I2))',err,error,*999)
508 IF(basis%NUMBER_OF_XI==3)
THEN 511 DO local_face_idx=1,basis%NUMBER_OF_LOCAL_FACES
514 & basis%LOCAL_FACE_XI_DIRECTION(local_face_idx),err,error,*999)
516 & basis%NUMBER_OF_NODES_IN_LOCAL_FACE(local_face_idx),err,error,*999)
518 & basis%NODE_NUMBERS_IN_LOCAL_FACE(:,local_face_idx),
'(" Nodes in local face :",4(X,I2))',
'(33X,4(X,I2))', &
520 DO nn=1,basis%NUMBER_OF_NODES_IN_LOCAL_FACE(local_face_idx)
523 & basis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE(1:,nn,local_face_idx),
'(" Derivatives in local face :",4(X,I2))', &
524 &
'(33X,4(X,I2))',err,error,*999)
530 DO local_line_idx=1,basis%NUMBER_OF_LOCAL_LINES
533 & basis%LOCAL_LINE_XI_DIRECTION(local_line_idx),err,error,*999)
535 & basis%NUMBER_OF_NODES_IN_LOCAL_LINE(local_line_idx),err,error,*999)
537 & basis%NODE_NUMBERS_IN_LOCAL_LINE(:,local_line_idx),
'(" Nodes in local line :",4(X,I2))',
'(33X,4(X,I2))', &
540 & basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE(:,local_line_idx),
'(" Derivatives in local line :",4(X,I2))', &
541 &
'(33X,4(X,I2))',err,error,*999)
542 IF(basis%NUMBER_OF_XI==2)
THEN 544 & basis%LOCAL_XI_NORMAL(local_line_idx),err,error,*999)
550 exits(
"BASIS_CREATE_FINISH")
552 999 errorsexits(
"BASIS_CREATE_FINISH",err,error)
577 INTEGER(INTG),
INTENT(IN) :: USER_NUMBER
579 INTEGER(INTG),
INTENT(OUT) :: ERR
589 enters(
"BASIS_CREATE_START",err,error,*999)
591 IF(
ASSOCIATED(basis))
THEN 592 CALL flagerror(
"Basis is already associated",err,error,*999)
596 IF(
ASSOCIATED(basis))
THEN 597 CALL flagerror(
"Basis number is already defined",err,error,*999)
601 IF(err/=0)
CALL flagerror(
"Could not allocate NEW_BASES",err,error,*999)
602 ALLOCATE(new_basis,stat=err)
603 IF(err/=0)
CALL flagerror(
"Could not allocate NEW_BASIS",err,error,*999)
613 new_basis%NUMBER_OF_SUB_BASES=0
614 NULLIFY(new_basis%SUB_BASES)
615 NULLIFY(new_basis%PARENT_BASIS)
617 new_basis%BASIS_FINISHED=.false.
618 new_basis%USER_NUMBER=user_number
619 new_basis%FAMILY_NUMBER=0
623 new_basis%NUMBER_OF_XI=3
624 ALLOCATE(new_basis%INTERPOLATION_XI(3),stat=err)
625 IF(err/=0)
CALL flagerror(
"Could not allocate basis interpolation xi",err,error,*999)
628 ALLOCATE(new_basis%COLLAPSED_XI(3),stat=err)
629 IF(err/=0)
CALL flagerror(
"Could not allocate basis collapsed xi",err,error,*999)
632 NULLIFY(new_basis%QUADRATURE%BASIS)
639 exits(
"BASIS_CREATE_START")
641 999
IF(
ASSOCIATED(new_basis))
CALL basis_destroy(new_basis,err,error,*998)
642 998
IF(
ASSOCIATED(new_bases))
DEALLOCATE(new_bases)
644 errorsexits(
"BASIS_CREATE_START",err,error)
657 INTEGER(INTG),
INTENT(IN) :: USER_NUMBER
658 INTEGER(INTG),
INTENT(OUT) :: ERR
662 enters(
"BASIS_DESTROY_NUMBER",err,error,*999)
666 exits(
"BASIS_DESTROY_NUMBER")
668 999 errorsexits(
"BASIS_DESTROY_NUMBER",err,error)
682 INTEGER(INTG),
INTENT(OUT) :: ERR
685 INTEGER(INTG) :: USER_NUMBER
687 enters(
"BASIS_DESTROY",err,error,*999)
689 IF(
ASSOCIATED(basis))
THEN 690 user_number=basis%USER_NUMBER
694 CALL flagerror(
"Basis is not associated.",err,error,*999)
697 exits(
"BASIS_DESTROY")
699 999 errorsexits(
"BASIS_DESTROY",err,error)
714 INTEGER(INTG),
INTENT(IN) :: ELEMENT_PARAMETER_INDEX
715 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIV_INDEX
716 REAL(DP),
INTENT(IN) :: XI(:)
717 INTEGER(INTG),
INTENT(OUT) :: ERR
720 REAL(DP) :: BASIS_EVALUATE_XI_DP
722 INTEGER(INTG) :: nn,nk
723 REAL(DP) :: XIL(size(xi,1)+1)
726 enters(
"BASIS_EVALUATE_XI_DP",err,error,*999)
728 basis_evaluate_xi_dp=0.0_dp
729 IF(
ASSOCIATED(basis))
THEN 730 IF(element_parameter_index>0.AND.element_parameter_index<=basis%NUMBER_OF_ELEMENT_PARAMETERS)
THEN 731 SELECT CASE(basis%TYPE)
733 nn=basis%ELEMENT_PARAMETER_INDEX_INV(1,element_parameter_index)
734 nk=basis%ELEMENT_PARAMETER_INDEX_INV(2,element_parameter_index)
739 xil(1:
SIZE(xi,1))=1.0_dp-xi
740 xil(
SIZE(xi,1)+1)=sum(xi)-(
SIZE(xi,1)-1.0_dp)
741 nn=basis%ELEMENT_PARAMETER_INDEX_INV(1,element_parameter_index)
749 local_error=
"Basis type "//
trim(
number_to_vstring(basis%TYPE,
"*",err,error))//
" is invalid or not implemented." 750 CALL flagerror(local_error,err,error,*999)
753 local_error=
"The specified element parameter index of "// &
755 &
" is invalid. The index must be > 0 and <= "// &
757 CALL flagerror(local_error,err,error,*999)
760 CALL flagerror(
"Basis is not associated.",err,error,*999)
763 exits(
"BASIS_EVALUATE_XI_DP")
765 999 errorsexits(
"BASIS_EVALUATE_XI_DP",err,error)
779 INTEGER(INTG),
INTENT(IN) :: order
780 INTEGER(INTG),
INTENT(IN) :: numCoords
781 INTEGER(INTG),
INTENT(OUT) :: numberGaussPoints
782 REAL(DP),
ALLOCATABLE,
INTENT(OUT) :: gaussPoints(:,:)
783 REAL(DP),
ALLOCATABLE,
INTENT(OUT) :: gaussWeights(:)
784 INTEGER(INTG),
INTENT(OUT) :: err
787 INTEGER(INTG) :: number_of_vertices,nj,ng,i,j,k,NUM_GAUSS_1,NUM_GAUSS_2,NUM_GAUSS_3,MAX_GAUSS
788 REAL(DP),
ALLOCATABLE :: XICOORDS(:,:),W(:,:),XI_MATRIX(:,:,:,:),XI(:)
791 enters(
"BASIS_GAUSS_POINTS_CALCULATE_DP",err,error,*999)
793 IF(
ASSOCIATED(basis))
THEN 795 SELECT CASE(numcoords)
805 max_gauss=order*order
810 max_gauss=order*order*order
812 local_error=
"Number of coordinates " &
814 CALL flagerror(local_error,err,error,*999)
817 ALLOCATE(w(max_gauss,numcoords),stat=err)
818 IF(err/=0)
CALL flagerror(
"Could not allocate weights",err,error,*999)
819 ALLOCATE(xi(numcoords),stat=err)
820 IF(err/=0)
CALL flagerror(
"Could not allocate gauss point coordinates",err,error,*999)
821 ALLOCATE(xicoords(max_gauss,numcoords),stat=err)
822 IF(err/=0)
CALL flagerror(
"Could not allocate gauss point coordinates",err,error,*999)
823 ALLOCATE(xi_matrix(max_gauss,max_gauss,max_gauss,numcoords),stat=err)
824 IF(err/=0)
CALL flagerror(
"Could not allocate XI matrix",err,error,*999)
826 SELECT CASE(basis%TYPE)
828 ALLOCATE(gausspoints(numcoords,max_gauss),stat=err)
829 IF(err/=0)
CALL flagerror(
"Could not allocate weights",err,error,*999)
830 ALLOCATE(gaussweights(max_gauss),stat=err)
831 IF(err/=0)
CALL flagerror(
"Could not allocate weights",err,error,*999)
833 CALL gauss_legendre(order,0.0_dp,1.0_dp,xicoords(1:order,nj),w(1:order,nj),err,error,*999)
841 xi_matrix(i,j,k,1)=xicoords(i,1)
842 xi_matrix(i,j,k,2)=xicoords(j,2)
843 xi_matrix(i,j,k,3)=xicoords(k,3)
844 xi(1:numcoords)=xi_matrix(i,j,k,1:numcoords)
845 ng=i+(j-1+(k-1)*num_gauss_2)*num_gauss_1
846 gaussweights(ng)=w(i,1)*w(j,2)*w(k,3)
847 gausspoints(1:numcoords,ng)=xi(1:numcoords)
848 numbergausspoints=numbergausspoints+1
853 IF(numcoords==1)
THEN 855 ELSEIF(numcoords==2)
THEN 860 ALLOCATE(gausspoints(number_of_vertices,max_gauss),stat=err)
861 IF(err/=0)
CALL flagerror(
"Could not allocate weights",err,error,*999)
862 ALLOCATE(gaussweights(max_gauss),stat=err)
863 IF(err/=0)
CALL flagerror(
"Could not allocate weights",err,error,*999)
865 CALL gauss_simplex(order,number_of_vertices,numbergausspoints,gausspoints,gaussweights,err,error,*999)
867 local_error=
"Basis type "// &
869 &
" is invalid or not implemented" 870 CALL flagerror(local_error,err,error,*999)
873 CALL flagerror(
"Basis is not associated",err,error,*999)
876 exits(
"BASIS_GAUSS_POINTS_CALCULATE_DP")
878 999 errorsexits(
"BASIS_GAUSS_POINTS_CALCULATE_DP",err,error)
892 INTEGER(INTG),
INTENT(IN) :: USER_NUMBER
893 INTEGER(INTG),
INTENT(IN) :: FAMILY_NUMBER
894 INTEGER(INTG),
INTENT(OUT) :: ERR
897 INTEGER(INTG) :: count,nb
901 enters(
"BASIS_FAMILY_DESTROY",err,error,*999)
905 IF(
ASSOCIATED(basis))
THEN 912 IF(basis%NUMBER_OF_SUB_BASES==0)
THEN 914 IF(
ASSOCIATED(basis%PARENT_BASIS))
THEN 916 NULLIFY(new_sub_bases)
917 IF(basis%PARENT_BASIS%NUMBER_OF_SUB_BASES>1)
THEN 919 ALLOCATE(new_sub_bases(basis%PARENT_BASIS%NUMBER_OF_SUB_BASES-1),stat=err)
920 IF(err/=0)
CALL flagerror(
"Could not allocate new sub-bases",err,error,*999)
922 DO nb=1,basis%PARENT_BASIS%NUMBER_OF_SUB_BASES
923 IF(basis%PARENT_BASIS%SUB_BASES(nb)%PTR%USER_NUMBER==basis%USER_NUMBER.AND. &
924 & basis%PARENT_BASIS%SUB_BASES(nb)%PTR%FAMILY_NUMBER/=basis%FAMILY_NUMBER)
THEN 926 new_sub_bases(count)%PTR=>basis%PARENT_BASIS%SUB_BASES(nb)%PTR
930 basis%PARENT_BASIS%NUMBER_OF_SUB_BASES=basis%PARENT_BASIS%NUMBER_OF_SUB_BASES-1
931 IF(
ASSOCIATED(basis%PARENT_BASIS%SUB_BASES))
DEALLOCATE(basis%PARENT_BASIS%SUB_BASES)
932 basis%PARENT_BASIS%SUB_BASES=>new_sub_bases
935 NULLIFY(new_sub_bases)
938 ALLOCATE(new_sub_bases(
basis_functions%NUMBER_BASIS_FUNCTIONS-1),stat=err)
939 IF(err/=0)
CALL flagerror(
"Could not allocate new sub-bases",err,error,*999)
958 DO WHILE(basis%NUMBER_OF_SUB_BASES>0)
959 CALL basis_family_destroy(basis%SUB_BASES(1)%PTR%USER_NUMBER,basis%SUB_BASES(1)%PTR%FAMILY_NUMBER,err,error,*999)
966 CALL flagerror(
"Basis user number does not exist",err,error,*999)
969 exits(
"BASIS_FAMILY_DESTROY")
971 999 errorsexits(
"BASIS_FAMILY_DESTROY",err,error)
984 INTEGER(INTG),
INTENT(IN) :: USER_NUMBER
985 INTEGER(INTG),
INTENT(IN) :: FAMILY_NUMBER
987 INTEGER(INTG),
INTENT(OUT) :: ERR
990 INTEGER(INTG) :: nb,nsb
993 enters(
"BASIS_FAMILY_NUMBER_FIND",err,error,*999)
997 DO WHILE(nb<=
basis_functions%NUMBER_BASIS_FUNCTIONS.AND..NOT.
ASSOCIATED(basis))
999 IF(family_number==0)
THEN 1004 DO WHILE(nsb<=
basis_functions%BASES(nb)%PTR%NUMBER_OF_SUB_BASES.AND..NOT.
ASSOCIATED(basis))
1006 IF(sub_basis%FAMILY_NUMBER==family_number)
THEN 1018 exits(
"BASIS_FAMILY_NUMBER_FIND")
1020 999 errorsexits(
"BASIS_FAMILY_NUMBER_FIND",err,error)
1033 INTEGER(INTG),
INTENT(OUT) :: ERR
1037 enters(
"BASIS_FINALISE",err,error,*999)
1039 IF(
ASSOCIATED(basis))
THEN 1040 IF(
ALLOCATED(basis%INTERPOLATION_XI))
DEALLOCATE(basis%INTERPOLATION_XI)
1041 IF(
ALLOCATED(basis%INTERPOLATION_TYPE))
DEALLOCATE(basis%INTERPOLATION_TYPE)
1042 IF(
ALLOCATED(basis%INTERPOLATION_ORDER))
DEALLOCATE(basis%INTERPOLATION_ORDER)
1043 IF(
ALLOCATED(basis%COLLAPSED_XI))
DEALLOCATE(basis%COLLAPSED_XI)
1044 IF(
ALLOCATED(basis%NODE_AT_COLLAPSE))
DEALLOCATE(basis%NODE_AT_COLLAPSE)
1046 IF(
ALLOCATED(basis%NUMBER_OF_NODES_XIC))
DEALLOCATE(basis%NUMBER_OF_NODES_XIC)
1047 IF(
ALLOCATED(basis%NUMBER_OF_DERIVATIVES))
DEALLOCATE(basis%NUMBER_OF_DERIVATIVES)
1048 IF(
ALLOCATED(basis%NODE_POSITION_INDEX))
DEALLOCATE(basis%NODE_POSITION_INDEX)
1049 IF(
ALLOCATED(basis%NODE_POSITION_INDEX_INV))
DEALLOCATE(basis%NODE_POSITION_INDEX_INV)
1050 IF(
ALLOCATED(basis%DERIVATIVE_ORDER_INDEX))
DEALLOCATE(basis%DERIVATIVE_ORDER_INDEX)
1051 IF(
ALLOCATED(basis%DERIVATIVE_ORDER_INDEX_INV))
DEALLOCATE(basis%DERIVATIVE_ORDER_INDEX_INV)
1052 IF(
ALLOCATED(basis%PARTIAL_DERIVATIVE_INDEX))
DEALLOCATE(basis%PARTIAL_DERIVATIVE_INDEX)
1053 IF(
ALLOCATED(basis%ELEMENT_PARAMETER_INDEX))
DEALLOCATE(basis%ELEMENT_PARAMETER_INDEX)
1054 IF(
ALLOCATED(basis%ELEMENT_PARAMETER_INDEX_INV))
DEALLOCATE(basis%ELEMENT_PARAMETER_INDEX_INV)
1055 IF(
ALLOCATED(basis%LOCAL_LINE_XI_DIRECTION))
DEALLOCATE(basis%LOCAL_LINE_XI_DIRECTION)
1056 IF(
ALLOCATED(basis%NUMBER_OF_NODES_IN_LOCAL_LINE))
DEALLOCATE(basis%NUMBER_OF_NODES_IN_LOCAL_LINE)
1057 IF(
ALLOCATED(basis%NODE_NUMBERS_IN_LOCAL_LINE))
DEALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_LINE)
1058 IF(
ALLOCATED(basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE))
DEALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE)
1059 IF(
ALLOCATED(basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE))
DEALLOCATE(basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE)
1060 IF(
ALLOCATED(basis%LOCAL_FACE_XI_DIRECTION))
DEALLOCATE(basis%LOCAL_FACE_XI_DIRECTION)
1061 IF(
ALLOCATED(basis%NUMBER_OF_NODES_IN_LOCAL_FACE))
DEALLOCATE(basis%NUMBER_OF_NODES_IN_LOCAL_FACE)
1062 IF(
ALLOCATED(basis%NODE_NUMBERS_IN_LOCAL_FACE))
DEALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_FACE)
1063 IF(
ALLOCATED(basis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE))
DEALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE)
1064 IF(
ALLOCATED(basis%ELEMENT_PARAMETERS_IN_LOCAL_FACE))
DEALLOCATE(basis%ELEMENT_PARAMETERS_IN_LOCAL_FACE)
1065 IF(
ALLOCATED(basis%LOCAL_XI_NORMAL))
DEALLOCATE(basis%LOCAL_XI_NORMAL)
1066 IF(
ASSOCIATED(basis%LINE_BASES))
DEALLOCATE(basis%LINE_BASES)
1067 IF(
ASSOCIATED(basis%FACE_BASES))
DEALLOCATE(basis%FACE_BASES)
1068 IF(
ASSOCIATED(basis%SUB_BASES))
DEALLOCATE(basis%SUB_BASES)
1072 exits(
"BASIS_FINALISE")
1074 999 errorsexits(
"BASIS_FINALISE",err,error)
1088 INTEGER(INTG),
INTENT(OUT) :: ERR
1092 enters(
"BASIS_INITIALISE",err,error,*999)
1094 IF(
ASSOCIATED(basis))
THEN 1096 basis%GLOBAL_NUMBER=0
1097 basis%FAMILY_NUMBER=0
1098 basis%BASIS_FINISHED=.false.
1099 basis%HERMITE=.false.
1101 basis%NUMBER_OF_XI=0
1102 basis%NUMBER_OF_XI_COORDINATES=0
1103 basis%DEGENERATE=.false.
1104 basis%NUMBER_OF_COLLAPSED_XI=0
1105 basis%NUMBER_OF_PARTIAL_DERIVATIVES=0
1106 basis%NUMBER_OF_NODES=0
1107 basis%NUMBER_OF_ELEMENT_PARAMETERS=0
1108 basis%MAXIMUM_NUMBER_OF_DERIVATIVES=0
1109 basis%NUMBER_OF_LOCAL_LINES=0
1110 basis%NUMBER_OF_LOCAL_FACES=0
1111 NULLIFY(basis%LINE_BASES)
1112 NULLIFY(basis%FACE_BASES)
1113 basis%NUMBER_OF_SUB_BASES=0
1114 NULLIFY(basis%SUB_BASES)
1115 NULLIFY(basis%PARENT_BASIS)
1117 CALL flagerror(
"Basis is not associated.",err,error,*999)
1120 exits(
"BASIS_INITIALISE")
1122 999 errorsexits(
"BASIS_INITIALISE",err,error)
1134 FUNCTION basis_interpolate_gauss_dp(BASIS,PARTIAL_DERIV_INDEX,QUADRATURE_SCHEME,GAUSS_POINT_NUMBER,ELEMENT_PARAMETERS,ERR,ERROR)
1138 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIV_INDEX
1139 INTEGER(INTG),
INTENT(IN) :: QUADRATURE_SCHEME
1140 INTEGER(INTG),
INTENT(IN) :: GAUSS_POINT_NUMBER
1141 REAL(DP),
INTENT(IN) :: ELEMENT_PARAMETERS(:)
1142 INTEGER(INTG),
INTENT(OUT) :: ERR
1145 REAL(DP) :: BASIS_INTERPOLATE_GAUSS_DP
1151 enters(
"BASIS_INTERPOLATE_GAUSS_DP",err,error,*999)
1153 basis_interpolate_gauss_dp=0.0_dp
1154 IF(
ASSOCIATED(basis))
THEN 1156 basis_quadrature_scheme=>basis%QUADRATURE%QUADRATURE_SCHEME_MAP(quadrature_scheme)%PTR
1157 IF(
ASSOCIATED(basis_quadrature_scheme))
THEN 1158 IF(gauss_point_number>0.AND.gauss_point_number<=basis_quadrature_scheme%NUMBER_OF_GAUSS)
THEN 1159 IF(partial_deriv_index>0.AND.partial_deriv_index<=basis%NUMBER_OF_PARTIAL_DERIVATIVES)
THEN 1160 DO ns=1,basis%NUMBER_OF_ELEMENT_PARAMETERS
1161 basis_interpolate_gauss_dp=basis_interpolate_gauss_dp+ &
1162 & basis_quadrature_scheme%GAUSS_BASIS_FNS(ns,partial_deriv_index,gauss_point_number)* &
1163 & element_parameters(ns)
1166 local_error=
"The partial derivative index of "//
trim(
number_to_vstring(partial_deriv_index,
"*",err,error))// &
1167 &
" is invalid. It must be between 1 and "// &
1169 CALL flagerror(local_error,err,error,*999)
1173 CALL flagerror(
"The quadrature scheme has not been created",err,error,*999)
1176 local_error=
"The quadrature scheme type number of "//
trim(
number_to_vstring(quadrature_scheme,
"*",err,error))// &
1177 &
" is invalid. It must be between 1 and "// &
1179 CALL flagerror(local_error,err,error,*999)
1182 CALL flagerror(
"Basis is not associated",err,error,*999)
1185 exits(
"BASIS_INTERPOLATE_GAUSS_DP")
1187 999 errorsexits(
"BASIS_INTERPOLATE_GAUSS_DP",err,error)
1198 & local_face_number,gauss_point_number,face_parameters,err,error)
1202 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIV_INDEX
1203 INTEGER(INTG),
INTENT(IN) :: QUADRATURE_SCHEME
1204 INTEGER(INTG),
INTENT(IN) :: LOCAL_FACE_NUMBER
1205 INTEGER(INTG),
INTENT(IN) :: GAUSS_POINT_NUMBER
1206 REAL(DP),
INTENT(IN) :: FACE_PARAMETERS(:)
1207 INTEGER(INTG),
INTENT(OUT) :: ERR
1210 REAL(DP) :: BASIS_INTERPOLATE_LOCAL_FACE_GAUSS_DP
1216 enters(
"BASIS_INTERPOLATE_LOCAL_FACE_GAUSS_DP",err,error,*999)
1218 basis_interpolate_local_face_gauss_dp=0.0_dp
1219 IF(
ASSOCIATED(basis))
THEN 1221 basis_quadrature_scheme=>basis%QUADRATURE%QUADRATURE_SCHEME_MAP(quadrature_scheme)%PTR
1222 IF(
ASSOCIATED(basis_quadrature_scheme))
THEN 1223 IF(basis%QUADRATURE%EVALUATE_FACE_GAUSS)
THEN 1224 IF(local_face_number>0.AND.local_face_number<=basis%NUMBER_OF_LOCAL_FACES)
THEN 1225 IF(gauss_point_number>0.AND.gauss_point_number<=basis_quadrature_scheme%NUMBER_OF_FACE_GAUSS(local_face_number))
THEN 1226 IF(partial_deriv_index>0.AND.partial_deriv_index<=basis%NUMBER_OF_PARTIAL_DERIVATIVES)
THEN 1227 DO ns=1,basis%NUMBER_OF_ELEMENT_PARAMETERS
1228 basis_interpolate_local_face_gauss_dp=basis_interpolate_local_face_gauss_dp+ &
1229 & basis_quadrature_scheme%FACE_GAUSS_BASIS_FNS(ns,partial_deriv_index,gauss_point_number,local_face_number)* &
1230 & face_parameters(ns)
1233 local_error=
"The partial derivative index of "//
trim(
number_to_vstring(partial_deriv_index,
"*",err,error))// &
1234 &
" is invalid. It must be between 1 and "// &
1236 CALL flagerror(local_error,err,error,*999)
1240 CALL flagerror(
"The local face number index is invalid.",err,error,*999)
1243 CALL flagerror(
"The face gauss interpolation scheme has not been created",err,error,*999)
1246 CALL flagerror(
"The quadrature scheme has not been created",err,error,*999)
1249 local_error=
"The quadrature scheme type number of "//
trim(
number_to_vstring(quadrature_scheme,
"*",err,error))// &
1250 &
" is invalid. It must be between 1 and "// &
1252 CALL flagerror(local_error,err,error,*999)
1255 CALL flagerror(
"Basis is not associated",err,error,*999)
1258 exits(
"BASIS_INTERPOLATE_LOCAL_FACE_GAUSS_DP")
1260 999 errorsexits(
"BASIS_INTERPOLATE_LOCAL_FACE_GAUSS_DP",err,error)
1275 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIV_INDEX
1276 REAL(DP),
INTENT(IN) :: XI(:)
1277 REAL(DP),
INTENT(IN) :: ELEMENT_PARAMETERS(:)
1278 INTEGER(INTG),
INTENT(OUT) :: ERR
1281 REAL(DP) :: BASIS_INTERPOLATE_XI_DP
1283 INTEGER(INTG) :: nn,nk,ns
1284 REAL(DP) :: XIL(size(xi,1)+1)
1287 enters(
"BASIS_INTERPOLATE_XI_DP",err,error,*999)
1289 basis_interpolate_xi_dp=0.0_dp
1290 IF(
ASSOCIATED(basis))
THEN 1291 SELECT CASE(basis%TYPE)
1294 DO nn=1,basis%NUMBER_OF_NODES
1295 DO nk=1,basis%NUMBER_OF_DERIVATIVES(nn)
1297 basis_interpolate_xi_dp=basis_interpolate_xi_dp+ &
1299 & element_parameters(ns)
1305 xil(1:
SIZE(xi,1))=1.0_dp-xi
1306 xil(
SIZE(xi,1)+1)=sum(xi)-(
SIZE(xi,1)-1.0_dp)
1308 DO nn=1,basis%NUMBER_OF_NODES
1310 basis_interpolate_xi_dp=basis_interpolate_xi_dp+ &
1312 & element_parameters(ns)
1319 local_error=
"Basis type "//
trim(
number_to_vstring(basis%TYPE,
"*",err,error))//
" is invalid or not implemented" 1320 CALL flagerror(local_error,err,error,*999)
1323 CALL flagerror(
"Basis is not associated",err,error,*999)
1326 exits(
"BASIS_INTERPOLATE_XI_DP")
1328 999 errorsexits(
"BASIS_INTERPOLATE_XI_DP",err,error)
1340 INTEGER(INTG),
INTENT(OUT) :: INTERPOLATION_XI(:)
1341 INTEGER(INTG),
INTENT(OUT) :: ERR
1346 enters(
"BASIS_INTERPOLATION_XI_GET",err,error,*999)
1348 IF(
ASSOCIATED(basis))
THEN 1349 IF(basis%BASIS_FINISHED)
THEN 1350 IF(
SIZE(interpolation_xi,1)>=
SIZE(basis%INTERPOLATION_XI,1))
THEN 1351 interpolation_xi=basis%INTERPOLATION_XI
1353 local_error=
"The size of INTERPOLATION_XI is too small. The supplied size is "// &
1356 CALL flagerror(local_error,err,error,*999)
1359 CALL flagerror(
"Basis has not been finished.",err,error,*999)
1362 CALL flagerror(
"Basis is not associated.",err,error,*999)
1365 exits(
"BASIS_INTERPOLATION_XI_GET")
1367 999 errorsexits(
"BASIS_INTERPOLATION_XI_GET",err,error)
1380 INTEGER(INTG),
INTENT(IN) :: USER_NUMBER
1381 INTEGER(INTG),
INTENT(IN) :: INTERPOLATION_XI(:)
1382 INTEGER(INTG),
INTENT(OUT) :: ERR
1387 enters(
"BASIS_INTERPOLATION_XI_SET_NUMBER",err,error,*999)
1392 exits(
"BASIS_INTERPOLATION_XI_SET_NUMBER")
1394 999 errorsexits(
"BASIS_INTERPOLATION_XI_SET_NUMBER",err,error)
1407 INTEGER(INTG),
INTENT(IN) :: INTERPOLATION_XI(:)
1408 INTEGER(INTG),
INTENT(OUT) :: ERR
1411 INTEGER(INTG) :: ni,LAST_INTERP
1414 enters(
"BASIS_INTERPOLATION_XI_SET_PTR",err,error,*999)
1416 IF(
ASSOCIATED(basis))
THEN 1417 IF(basis%BASIS_FINISHED)
THEN 1418 CALL flagerror(
"Basis has been finished",err,error,*999)
1420 IF(
SIZE(interpolation_xi,1)==basis%NUMBER_OF_XI)
THEN 1422 SELECT CASE(basis%TYPE)
1424 DO ni=1,basis%NUMBER_OF_XI
1425 SELECT CASE(interpolation_xi(ni))
1430 local_error=
"Interpolation xi value "//
trim(
number_to_vstring(interpolation_xi(ni),
"*",err,error))// &
1431 &
" for xi direction "//
trim(
number_to_vstring(ni,
"*",err,error))//
" is invalid for a Lagrange-Hermite TP basis." 1432 CALL flagerror(local_error,err,error,*999)
1436 last_interp=interpolation_xi(1)
1437 DO ni=1,basis%NUMBER_OF_XI
1438 SELECT CASE(interpolation_xi(ni))
1440 IF(interpolation_xi(ni)/=last_interp)
THEN 1441 CALL flagerror(
"The interpolation xi value must be the same for all xi directions for a simplex basis.", &
1445 local_error=
"Interpolation xi value "//
trim(
number_to_vstring(interpolation_xi(ni),
"*",err,error))// &
1446 &
" for xi direction "//
trim(
number_to_vstring(ni,
"*",err,error))//
" is invalid for a simplex basis." 1447 CALL flagerror(local_error,err,error,*999)
1451 CALL flagerror(
"Invalid basis type or not implemented",err,error,*999)
1454 basis%INTERPOLATION_XI(1:basis%NUMBER_OF_XI)=interpolation_xi(1:basis%NUMBER_OF_XI)
1456 local_error=
"The size of the interpolation xi array ("// &
1457 &
trim(
number_to_vstring(
SIZE(interpolation_xi,1),
"*",err,error))//
") does not match the number of xi directions ("// &
1460 CALL flagerror(local_error,err,error,*999)
1464 CALL flagerror(
"Basis is not associated.",err,error,*999)
1467 exits(
"BASIS_INTERPOLATION_XI_SET_PTR")
1469 999 errorsexits(
"BASIS_INTERPOLATION_XI_SET_PTR",err,error)
1483 INTEGER(INTG),
INTENT(OUT) :: err
1486 INTEGER(INTG) :: maximumNumberOfNodes,numberOfDerivatives,xiIdx,xiIdx1,xiIdx2,xiIdx3,derivativeIdx,localNode,localLineNodeIdx, &
1487 & elementParameter,oldNumberOfDerivatives,position(4),collapsedPosition(3),maximumNodeExtent(3),collapsedXi(3), &
1488 & numberOfNodes,numberOfLocalLines,nodeCount,specialNodeCount,nodesInLine(4),numberOfLocalFaces,localFaceIdx, &
1489 & localNodeIdx,localNodeIdx1,localNodeIdx2,localNodeIdx3,directionIdx,localFaceDerivative,localNodeCount, &
1490 & localLineParameter,localFaceParameter
1491 LOGICAL,
ALLOCATABLE :: nodeAtCollapse(:)
1492 LOGICAL :: atCollapse,collapsedFace,firstCollapsedPosition
1494 enters(
"Basis_LHTPBasisCreate",err,error,*999)
1496 IF(
ASSOCIATED(basis))
THEN 1498 basis%NUMBER_OF_XI_COORDINATES=basis%NUMBER_OF_XI
1499 basis%NUMBER_OF_PARTIAL_DERIVATIVES=basis%NUMBER_OF_XI_COORDINATES**2+2
1500 ALLOCATE(basis%INTERPOLATION_TYPE(basis%NUMBER_OF_XI_COORDINATES),stat=err)
1501 IF(err/=0)
CALL flagerror(
"Could not allocate basis interpolation type.",err,error,*999)
1502 ALLOCATE(basis%INTERPOLATION_ORDER(basis%NUMBER_OF_XI_COORDINATES),stat=err)
1503 IF(err/=0)
CALL flagerror(
"Could not allocate basis interpolation order.",err,error,*999)
1504 ALLOCATE(basis%NUMBER_OF_NODES_XIC(basis%NUMBER_OF_XI_COORDINATES),stat=err)
1505 IF(err/=0)
CALL flagerror(
"Could not allocate basis number of nodes xic.",err,error,*999)
1507 maximumnumberofnodes=0
1508 basis%degenerate=.false.
1509 basis%NUMBER_OF_COLLAPSED_XI=0
1510 DO xiidx=1,basis%NUMBER_OF_XI
1512 SELECT CASE(basis%INTERPOLATION_XI(xiidx))
1516 basis%NUMBER_OF_NODES_XIC(xiidx)=2
1520 basis%NUMBER_OF_NODES_XIC(xiidx)=3
1524 basis%NUMBER_OF_NODES_XIC(xiidx)=4
1528 basis%NUMBER_OF_NODES_XIC(xiidx)=2
1532 basis%NUMBER_OF_NODES_XIC(xiidx)=2
1536 basis%NUMBER_OF_NODES_XIC(xiidx)=2
1538 CALL flagerror(
"Invalid interpolation type",err,error,*999)
1541 basis%NUMBER_OF_COLLAPSED_XI=basis%NUMBER_OF_COLLAPSED_XI+1
1542 collapsedxi(basis%NUMBER_OF_COLLAPSED_XI)=xiidx
1543 basis%degenerate=.true.
1545 numberofnodes=numberofnodes*basis%NUMBER_OF_NODES_XIC(xiidx)
1546 IF(basis%NUMBER_OF_NODES_XIC(xiidx)>maximumnumberofnodes) maximumnumberofnodes=basis%NUMBER_OF_NODES_XIC(xiidx)
1549 IF(basis%degenerate)
THEN 1551 ALLOCATE(nodeatcollapse(numberofnodes),stat=err)
1552 IF(err/=0)
CALL flagerror(
"Could not allocate at collapse",err,error,*999)
1554 basis%NUMBER_OF_NODES=0
1556 DO localnodeidx=1,numberofnodes
1558 DO xiidx=1,basis%NUMBER_OF_XI
1560 & basis%COLLAPSED_XI(xiidx)==
basis_collapsed_at_xi1.AND.position(xiidx)==basis%NUMBER_OF_NODES_XIC(xiidx)) &
1563 firstcollapsedposition=all(position(collapsedxi(1:basis%NUMBER_OF_COLLAPSED_XI))==1)
1568 IF(firstcollapsedposition)
THEN 1569 basis%NUMBER_OF_NODES=basis%NUMBER_OF_NODES+1
1570 nodeatcollapse(basis%NUMBER_OF_NODES)=.true.
1573 basis%NUMBER_OF_NODES=basis%NUMBER_OF_NODES+1
1574 nodeatcollapse(basis%NUMBER_OF_NODES)=.false.
1576 position(1)=position(1)+1
1577 DO xiidx=1,basis%NUMBER_OF_XI
1578 IF(position(xiidx)>basis%NUMBER_OF_NODES_XIC(xiidx))
THEN 1580 position(xiidx+1)=position(xiidx+1)+1
1584 CALL move_alloc(nodeatcollapse,basis%NODE_AT_COLLAPSE)
1586 basis%NUMBER_OF_NODES=numberofnodes
1587 ALLOCATE(basis%NODE_AT_COLLAPSE(basis%NUMBER_OF_NODES),stat=err)
1588 IF(err/=0)
CALL flagerror(
"Could not allocate basis node at collapse.",err,error,*999)
1589 basis%NODE_AT_COLLAPSE=.false.
1593 ALLOCATE(basis%NODE_POSITION_INDEX(basis%NUMBER_OF_NODES,basis%NUMBER_OF_XI_COORDINATES),stat=err)
1594 IF(err/=0)
CALL flagerror(
"Could not allocate basis node position index.",err,error,*999)
1595 SELECT CASE(basis%NUMBER_OF_XI_COORDINATES)
1597 ALLOCATE(basis%NODE_POSITION_INDEX_INV(maximumnumberofnodes,1,1,1),stat=err)
1599 ALLOCATE(basis%NODE_POSITION_INDEX_INV(maximumnumberofnodes,maximumnumberofnodes,1,1),stat=err)
1601 ALLOCATE(basis%NODE_POSITION_INDEX_INV(maximumnumberofnodes,maximumnumberofnodes,maximumnumberofnodes,1),stat=err)
1603 CALL flagerror(
"Invalid number of xi coordinates.",err,error,*999)
1605 IF(err/=0)
CALL flagerror(
"Could not allocate node position index inverse.",err,error,*999)
1606 basis%NODE_POSITION_INDEX_INV=0
1612 firstcollapsedposition=.true.
1613 DO localnodeidx1=1,numberofnodes
1615 IF(basis%degenerate)
THEN 1616 DO xiidx=1,basis%NUMBER_OF_XI
1618 & basis%COLLAPSED_XI(xiidx)==
basis_collapsed_at_xi1.AND.position(xiidx)==basis%NUMBER_OF_NODES_XIC(xiidx)) &
1621 firstcollapsedposition=all(position(collapsedxi(1:basis%NUMBER_OF_COLLAPSED_XI))==1)
1627 IF(firstcollapsedposition)
THEN 1628 localnode=localnode+1
1629 basis%NODE_POSITION_INDEX(localnode,1:basis%NUMBER_OF_XI)=position(1:basis%NUMBER_OF_XI)
1630 basis%NODE_POSITION_INDEX_INV(position(1),position(2),position(3),1)=localnode
1633 collapsedposition(1:basis%NUMBER_OF_XI)=position(1:basis%NUMBER_OF_XI)
1634 collapsedposition(collapsedxi(1:basis%NUMBER_OF_COLLAPSED_XI))=1
1635 basis%NODE_POSITION_INDEX_INV(position(1),position(2),position(3),1)= &
1636 & basis%NODE_POSITION_INDEX_INV(collapsedposition(1),collapsedposition(2),collapsedposition(3),1)
1639 localnode=localnode+1
1640 basis%NODE_POSITION_INDEX(localnode,1:basis%NUMBER_OF_XI)=position(1:basis%NUMBER_OF_XI)
1641 basis%NODE_POSITION_INDEX_INV(position(1),position(2),position(3),1)=localnode
1643 position(1)=position(1)+1
1644 DO xiidx=1,basis%NUMBER_OF_XI
1645 IF(position(xiidx)>basis%NUMBER_OF_NODES_XIC(xiidx))
THEN 1647 position(xiidx+1)=position(xiidx+1)+1
1652 basis%MAXIMUM_NUMBER_OF_DERIVATIVES=-1
1653 basis%NUMBER_OF_ELEMENT_PARAMETERS=0
1654 DO localnodeidx=1,basis%NUMBER_OF_NODES
1655 numberofderivatives=1
1656 DO xiidx=1,basis%NUMBER_OF_XI
1657 IF((.NOT.basis%NODE_AT_COLLAPSE(localnodeidx).OR.basis%COLLAPSED_XI(xiidx)==
basis_not_collapsed).AND. &
1660 & (basis%NODE_POSITION_INDEX(localnodeidx,xiidx)==1.AND. &
1662 & (basis%NODE_POSITION_INDEX(localnodeidx,xiidx)==2.AND. &
1665 numberofderivatives=numberofderivatives*2
1668 basis%NUMBER_OF_ELEMENT_PARAMETERS=basis%NUMBER_OF_ELEMENT_PARAMETERS+numberofderivatives
1669 IF(numberofderivatives>basis%MAXIMUM_NUMBER_OF_DERIVATIVES) basis%MAXIMUM_NUMBER_OF_DERIVATIVES=numberofderivatives
1672 ALLOCATE(basis%NUMBER_OF_DERIVATIVES(basis%NUMBER_OF_NODES),stat=err)
1673 IF(err/=0)
CALL flagerror(
"Could not allocate number of derivatives.",err,error,*999)
1674 ALLOCATE(basis%DERIVATIVE_ORDER_INDEX(basis%MAXIMUM_NUMBER_OF_DERIVATIVES,basis%NUMBER_OF_NODES, &
1675 & basis%NUMBER_OF_XI),stat=err)
1676 IF(err/=0)
CALL flagerror(
"Could not allocate derivative order index.",err,error,*999)
1678 & basis%NUMBER_OF_NODES),stat=err)
1679 IF(err/=0)
CALL flagerror(
"Could not allocate derivative order index inverse.",err,error,*999)
1680 ALLOCATE(basis%PARTIAL_DERIVATIVE_INDEX(basis%MAXIMUM_NUMBER_OF_DERIVATIVES,basis%NUMBER_OF_NODES),stat=err)
1681 IF(err/=0)
CALL flagerror(
"Could not allocate partial derivative index.",err,error,*999)
1682 ALLOCATE(basis%ELEMENT_PARAMETER_INDEX(basis%MAXIMUM_NUMBER_OF_DERIVATIVES,basis%NUMBER_OF_NODES),stat=err)
1683 IF(err/=0)
CALL flagerror(
"Could not allocate element parameter index.",err,error,*999)
1684 ALLOCATE(basis%ELEMENT_PARAMETER_INDEX_INV(2,basis%NUMBER_OF_ELEMENT_PARAMETERS),stat=err)
1685 IF(err/=0)
CALL flagerror(
"Could not allocate element parameter index inverse.",err,error,*999)
1688 basis%DERIVATIVE_ORDER_INDEX=0
1689 basis%DERIVATIVE_ORDER_INDEX_INV=0
1690 DO localnodeidx=1,basis%NUMBER_OF_NODES
1691 basis%NUMBER_OF_DERIVATIVES(localnodeidx)=1
1692 DO xiidx1=1,basis%NUMBER_OF_XI
1693 IF((.NOT.basis%NODE_AT_COLLAPSE(localnodeidx).OR.basis%COLLAPSED_XI(xiidx1)==
basis_not_collapsed).AND. &
1696 & (basis%NODE_POSITION_INDEX(localnodeidx,xiidx1)==1.AND. &
1698 & (basis%NODE_POSITION_INDEX(localnodeidx,xiidx1)==2.AND. &
1700 oldnumberofderivatives=basis%NUMBER_OF_DERIVATIVES(localnodeidx)
1701 basis%NUMBER_OF_DERIVATIVES(localnodeidx)=basis%NUMBER_OF_DERIVATIVES(localnodeidx)*2
1702 DO derivativeidx=1,oldnumberofderivatives
1703 basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,xiidx1)=
no_part_deriv 1704 basis%DERIVATIVE_ORDER_INDEX(oldnumberofderivatives+derivativeidx,localnodeidx,xiidx1)=
first_part_deriv 1705 DO xiidx2=1,xiidx1-1
1706 basis%DERIVATIVE_ORDER_INDEX(oldnumberofderivatives+derivativeidx,localnodeidx,xiidx2)= &
1707 & basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,xiidx2)
1711 DO derivativeidx=1,basis%NUMBER_OF_DERIVATIVES(localnodeidx)
1712 basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,xiidx1)=
no_part_deriv 1717 DO derivativeidx=1,basis%NUMBER_OF_DERIVATIVES(localnodeidx)
1718 elementparameter=elementparameter+1
1719 basis%ELEMENT_PARAMETER_INDEX(derivativeidx,localnodeidx)=elementparameter
1720 basis%ELEMENT_PARAMETER_INDEX_INV(1,elementparameter)=localnodeidx
1721 basis%ELEMENT_PARAMETER_INDEX_INV(2,elementparameter)=derivativeidx
1722 SELECT CASE(basis%NUMBER_OF_XI)
1724 basis%DERIVATIVE_ORDER_INDEX_INV(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,1),1,1,localnodeidx)= &
1726 SELECT CASE(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,1))
1728 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
no_part_deriv 1730 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
part_deriv_s1 1732 CALL flagerror(
"Invalid derivative order index.",err,error,*999)
1735 basis%DERIVATIVE_ORDER_INDEX_INV(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,1), &
1736 & basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,2),1,localnodeidx)=derivativeidx
1737 SELECT CASE(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,1))
1739 SELECT CASE(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,2))
1741 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
no_part_deriv 1743 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
part_deriv_s2 1745 CALL flagerror(
"Invalid derivative order index.",err,error,*999)
1748 SELECT CASE(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,2))
1750 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
part_deriv_s1 1752 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
part_deriv_s1_s2 1754 CALL flagerror(
"Invalid derivative order index.",err,error,*999)
1757 CALL flagerror(
"Invalid derivative order index.",err,error,*999)
1760 basis%DERIVATIVE_ORDER_INDEX_INV(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,1), &
1761 & basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,2), &
1762 & basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,3),localnodeidx)=derivativeidx
1763 SELECT CASE(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,1))
1765 SELECT CASE(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,2))
1767 SELECT CASE(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,3))
1769 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
no_part_deriv 1771 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
part_deriv_s3 1773 CALL flagerror(
"Invalid derivative order index.",err,error,*999)
1776 SELECT CASE(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,3))
1778 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
part_deriv_s2 1780 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
part_deriv_s2_s3 1782 CALL flagerror(
"Invalid derivative order index.",err,error,*999)
1785 CALL flagerror(
"Invalid derivative order index.",err,error,*999)
1788 SELECT CASE(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,2))
1790 SELECT CASE(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,3))
1792 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
part_deriv_s1 1794 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
part_deriv_s1_s3 1796 CALL flagerror(
"Invalid derivative order index.",err,error,*999)
1799 SELECT CASE(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx,3))
1801 basis%PARTIAL_DERIVATIVE_INDEX(derivativeidx,localnodeidx)=
part_deriv_s1_s2 1805 CALL flagerror(
"Invalid derivative order index.",err,error,*999)
1808 CALL flagerror(
"Invalid derivative order index.",err,error,*999)
1811 CALL flagerror(
"Invalid derivative order index.",err,error,*999)
1814 CALL flagerror(
"Invalid number of xi direcions.",err,error,*999)
1820 SELECT CASE(basis%NUMBER_OF_XI)
1822 numberoflocallines=1
1823 basis%NUMBER_OF_LOCAL_LINES=1
1824 ALLOCATE(basis%NUMBER_OF_NODES_IN_LOCAL_LINE(numberoflocallines),stat=err)
1825 IF(err/=0)
CALL flagerror(
"Could not allocate number of nodes in local line.",err,error,*999)
1826 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(1)=basis%NUMBER_OF_NODES_XIC(1)
1827 ALLOCATE(basis%LOCAL_LINE_XI_DIRECTION(numberoflocallines),stat=err)
1828 IF(err/=0)
CALL flagerror(
"Could not allocate local line xi direction.",err,error,*999)
1829 basis%LOCAL_LINE_XI_DIRECTION(1)=1
1830 ALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_LINE(basis%NUMBER_OF_NODES_XIC(1),numberoflocallines),stat=err)
1831 IF(err/=0)
CALL flagerror(
"Could not allocate node numbers in local line.",err,error,*999)
1832 ALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE(basis%NUMBER_OF_NODES_XIC(1),numberoflocallines),stat=err)
1833 IF(err/=0)
CALL flagerror(
"Could not allocate derivative numbers in local line.",err,error,*999)
1835 ALLOCATE(basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE(basis%NUMBER_OF_NODES_XIC(1)**2,numberoflocallines),stat=err)
1836 IF(err/=0)
CALL flagerror(
"Could not allocate element parameters in local line.",err,error,*999)
1837 basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE=1
1838 locallineparameter=0
1839 DO localnodeidx2=1,basis%NUMBER_OF_NODES_XIC(1)
1840 DO localnodeidx1=1,basis%NUMBER_OF_NODES
1841 IF(basis%NODE_POSITION_INDEX(localnodeidx1,1)==localnodeidx2)
THEN 1842 basis%NODE_NUMBERS_IN_LOCAL_LINE(localnodeidx2,1)=localnodeidx1
1843 DO derivativeidx=1,basis%NUMBER_OF_DERIVATIVES(localnodeidx2)
1844 locallineparameter=locallineparameter+1
1845 basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE(locallineparameter,1)=basis%ELEMENT_PARAMETER_INDEX( &
1846 & derivativeidx,localnodeidx1)
1847 IF(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnodeidx2,1)==
first_part_deriv)
THEN 1848 basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE(localnodeidx2,1)=derivativeidx
1858 maximumnodeextent(1)=maxval(basis%NODE_POSITION_INDEX(:,1))
1859 maximumnodeextent(2)=maxval(basis%NODE_POSITION_INDEX(:,2))
1861 numberoflocallines=4-basis%NUMBER_OF_COLLAPSED_XI
1862 ALLOCATE(basis%NUMBER_OF_NODES_IN_LOCAL_LINE(numberoflocallines),stat=err)
1863 IF(err/=0)
CALL flagerror(
"Could not allocate number of nodes in local line.",err,error,*999)
1864 basis%NUMBER_OF_NODES_IN_LOCAL_LINE=0
1865 ALLOCATE(basis%LOCAL_LINE_XI_DIRECTION(numberoflocallines),stat=err)
1866 IF(err/=0)
CALL flagerror(
"Could not allocate local line xi direction.",err,error,*999)
1867 ALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC),numberoflocallines),stat=err)
1868 IF(err/=0)
CALL flagerror(
"Could not allocate node numbers in local line",err,error,*999)
1869 basis%NODE_NUMBERS_IN_LOCAL_LINE=0
1870 ALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC),numberoflocallines),stat=err)
1871 IF(err/=0)
CALL flagerror(
"Could not allocate derivative numbers in local line.",err,error,*999)
1873 ALLOCATE(basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC)**2,numberoflocallines),stat=err)
1874 IF(err/=0)
CALL flagerror(
"Could not allocate element parameters in local line.",err,error,*999)
1875 basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE=1
1876 ALLOCATE(basis%LOCAL_XI_NORMAL(numberoflocallines),stat=err)
1877 IF(err/=0)
CALL flagerror(
"Could not allocate local xi normal.",err,error,*999)
1879 basis%NUMBER_OF_LOCAL_LINES=0
1884 DO localnodeidx2=1,maximumnodeextent(xiidx2),maximumnodeextent(xiidx2)-1
1888 DO localnodeidx1=1,basis%NUMBER_OF_NODES
1894 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2.OR. &
1895 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx1)==1)
THEN 1896 nodecount=nodecount+1
1897 nodesinline(nodecount)=localnodeidx1
1900 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2)
THEN 1901 nodecount=nodecount+1
1902 nodesinline(nodecount)=localnodeidx1
1903 ELSE IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx1)==maximumnodeextent(xiidx1))
THEN 1905 specialnodecount=specialnodecount+1
1906 nodesinline(maximumnodeextent(xiidx1))=localnodeidx1
1908 nodecount=nodecount+1
1909 nodesinline(nodecount)=localnodeidx1
1915 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2)
THEN 1916 nodecount=nodecount+1
1917 nodesinline(nodecount)=localnodeidx1
1922 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2)
THEN 1923 nodecount=nodecount+1
1924 nodesinline(nodecount)=localnodeidx1
1928 IF((nodecount+specialnodecount)>1)
THEN 1929 basis%NUMBER_OF_LOCAL_LINES=basis%NUMBER_OF_LOCAL_LINES+1
1930 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(basis%NUMBER_OF_LOCAL_LINES)=nodecount+specialnodecount
1931 basis%NODE_NUMBERS_IN_LOCAL_LINE(1:basis%NUMBER_OF_NODES_IN_LOCAL_LINE(basis%NUMBER_OF_LOCAL_LINES), &
1932 & basis%NUMBER_OF_LOCAL_LINES)=nodesinline(1:basis%NUMBER_OF_NODES_IN_LOCAL_LINE(basis%NUMBER_OF_LOCAL_LINES))
1933 locallineparameter=0
1934 DO locallinenodeidx=1,basis%NUMBER_OF_NODES_IN_LOCAL_LINE(basis%NUMBER_OF_LOCAL_LINES)
1935 DO derivativeidx=1,basis%NUMBER_OF_DERIVATIVES(basis%NODE_NUMBERS_IN_LOCAL_LINE(locallinenodeidx, &
1936 & basis%NUMBER_OF_LOCAL_LINES))
1937 IF(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,basis%NODE_NUMBERS_IN_LOCAL_LINE(locallinenodeidx, &
1939 locallineparameter=locallineparameter+1
1940 basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE(locallineparameter,basis%NUMBER_OF_LOCAL_LINES)= &
1941 & basis%ELEMENT_PARAMETER_INDEX(derivativeidx,basis%NODE_NUMBERS_IN_LOCAL_LINE(locallinenodeidx, &
1942 & basis%NUMBER_OF_LOCAL_LINES))
1943 IF(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,basis%NODE_NUMBERS_IN_LOCAL_LINE(locallinenodeidx, &
1945 basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE(locallinenodeidx,basis%NUMBER_OF_LOCAL_LINES)=derivativeidx
1950 basis%LOCAL_LINE_XI_DIRECTION(basis%NUMBER_OF_LOCAL_LINES)=xiidx1
1951 IF(localnodeidx2==1)
THEN 1952 basis%LOCAL_XI_NORMAL(basis%NUMBER_OF_LOCAL_LINES)=-xiidx2
1954 basis%LOCAL_XI_NORMAL(basis%NUMBER_OF_LOCAL_LINES)=xiidx2
1961 maximumnodeextent(1)=maxval(basis%NODE_POSITION_INDEX(:,1))
1962 maximumnodeextent(2)=maxval(basis%NODE_POSITION_INDEX(:,2))
1963 maximumnodeextent(3)=maxval(basis%NODE_POSITION_INDEX(:,3))
1965 IF(basis%NUMBER_OF_COLLAPSED_XI==1)
THEN 1966 numberoflocallines=9
1967 numberoflocalfaces=5
1968 ELSE IF(basis%NUMBER_OF_COLLAPSED_XI==2)
THEN 1969 numberoflocallines=8
1970 numberoflocalfaces=5
1972 numberoflocallines=12
1973 numberoflocalfaces=6
1975 basis%NUMBER_OF_LOCAL_FACES=numberoflocalfaces
1977 ALLOCATE(basis%NUMBER_OF_NODES_IN_LOCAL_LINE(numberoflocallines),stat=err)
1978 IF(err/=0)
CALL flagerror(
"Could not allocate number of nodes in local line.",err,error,*999)
1979 basis%NUMBER_OF_NODES_IN_LOCAL_LINE=0
1981 ALLOCATE(basis%NUMBER_OF_NODES_IN_LOCAL_FACE(numberoflocalfaces),stat=err)
1982 IF(err/=0)
CALL flagerror(
"Could not allocate number of nodes in local face.",err,error,*999)
1983 basis%NUMBER_OF_NODES_IN_LOCAL_FACE=0
1985 ALLOCATE(basis%LOCAL_LINE_XI_DIRECTION(numberoflocallines),stat=err)
1986 IF(err/=0)
CALL flagerror(
"Could not allocate local line xi direction.",err,error,*999)
1988 ALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC),numberoflocallines),stat=err)
1989 IF(err/=0)
CALL flagerror(
"Could not allocate node numbers in local line.",err,error,*999)
1990 basis%NODE_NUMBERS_IN_LOCAL_LINE=0
1992 ALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC),numberoflocallines),stat=err)
1993 IF(err/=0)
CALL flagerror(
"Could not allocate derivative numbers in local line.",err,error,*999)
1996 ALLOCATE(basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC)**2,numberoflocallines),stat=err)
1997 IF(err/=0)
CALL flagerror(
"Could not allocate element parameters in local line.",err,error,*999)
1998 basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE=1
2000 ALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE(0:basis%MAXIMUM_NUMBER_OF_DERIVATIVES, &
2001 & maxval(basis%NUMBER_OF_NODES_XIC)**2,numberoflocalfaces),stat=err)
2002 IF(err/=0)
CALL flagerror(
"Could not allocate derivative numbers in local face.",err,error,*999)
2004 basis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE(0,:,:)=1
2006 ALLOCATE(basis%ELEMENT_PARAMETERS_IN_LOCAL_FACE(maxval(basis%NUMBER_OF_NODES_XIC)**2* &
2007 & basis%MAXIMUM_NUMBER_OF_DERIVATIVES,numberoflocalfaces),stat=err)
2008 IF(err/=0)
CALL flagerror(
"Could not allocate element parameters in local face.",err,error,*999)
2009 basis%ELEMENT_PARAMETERS_IN_LOCAL_FACE=1
2011 ALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_FACE(max(maximumnodeextent(2)*maximumnodeextent(3), &
2012 & maximumnodeextent(3)*maximumnodeextent(1),maximumnodeextent(2)*maximumnodeextent(1)), &
2013 & numberoflocalfaces),stat=err)
2014 IF(err/=0)
CALL flagerror(
"Could not allocate node numbers in local face.",err,error,*999)
2015 basis%NODE_NUMBERS_IN_LOCAL_FACE=0
2017 ALLOCATE(basis%LOCAL_XI_NORMAL(numberoflocalfaces),stat=err)
2018 IF(err/=0)
CALL flagerror(
"Could not allocate local xi normal.",err,error,*999)
2020 ALLOCATE(basis%LOCAL_FACE_XI_DIRECTION(numberoflocalfaces),stat=err)
2021 IF(err/=0)
CALL flagerror(
"Could not allocate local face xi direction.",err,error,*999)
2024 basis%NUMBER_OF_LOCAL_LINES=0
2029 DO localnodeidx3=1,maximumnodeextent(xiidx3),maximumnodeextent(xiidx3)-1
2030 DO localnodeidx2=1,maximumnodeextent(xiidx2),maximumnodeextent(xiidx2)-1
2035 DO localnodeidx1=1,basis%NUMBER_OF_NODES
2041 IF((basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2.OR. &
2042 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx1)==1).AND. &
2043 & (basis%NODE_POSITION_INDEX(localnodeidx1,xiidx3)==localnodeidx3.OR. &
2044 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx1)==1))
THEN 2045 nodecount=nodecount+1
2046 nodesinline(nodecount)=localnodeidx1
2049 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2.AND. &
2050 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx3)==localnodeidx3)
THEN 2051 nodecount=nodecount+1
2052 nodesinline(nodecount)=localnodeidx1
2053 ELSE IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx1)==maximumnodeextent(xiidx1))
THEN 2055 specialnodecount=specialnodecount+1
2056 nodesinline(maximumnodeextent(xiidx1))=localnodeidx1
2058 nodecount=nodecount+1
2059 nodesinline(nodecount)=localnodeidx1
2067 IF((basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2.OR. &
2068 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx1)==1).AND. &
2069 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx3)==localnodeidx3)
THEN 2070 nodecount=nodecount+1
2071 nodesinline(nodecount)=localnodeidx1
2074 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2.AND. &
2075 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx3)==localnodeidx3)
THEN 2076 nodecount=nodecount+1
2077 nodesinline(nodecount)=localnodeidx1
2078 ELSE IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx1)==maximumnodeextent(xiidx1).AND. &
2079 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx3)==localnodeidx3)
THEN 2080 IF(xiidx1<xiidx2)
THEN 2081 specialnodecount=specialnodecount+1
2082 nodesinline(maximumnodeextent(xiidx1))=localnodeidx1
2084 nodecount=nodecount+1
2085 nodesinline(nodecount)=localnodeidx1
2090 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2.AND. &
2091 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx3)==localnodeidx3)
THEN 2092 nodecount=nodecount+1
2093 nodesinline(nodecount)=localnodeidx1
2099 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2.AND. &
2100 & (basis%NODE_POSITION_INDEX(localnodeidx1,xiidx3)==localnodeidx3.OR. &
2101 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx1)==1))
THEN 2102 nodecount=nodecount+1
2103 nodesinline(nodecount)=localnodeidx1
2106 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2.AND. &
2107 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx3)==localnodeidx3)
THEN 2108 nodecount=nodecount+1
2109 nodesinline(nodecount)=localnodeidx1
2110 ELSE IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx1)==maximumnodeextent(xiidx1).AND. &
2111 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2)
THEN 2112 IF(xiidx1<xiidx3)
THEN 2113 specialnodecount=specialnodecount+1
2114 nodesinline(maximumnodeextent(xiidx1))=localnodeidx1
2116 nodecount=nodecount+1
2117 nodesinline(nodecount)=localnodeidx1
2122 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2.AND. &
2123 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx3)==localnodeidx3)
THEN 2124 nodecount=nodecount+1
2125 nodesinline(nodecount)=localnodeidx1
2130 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2.AND. &
2131 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx3)==localnodeidx3)
THEN 2132 nodecount=nodecount+1
2133 nodesinline(nodecount)=localnodeidx1
2139 IF(basis%NODE_POSITION_INDEX(localnodeidx1,xiidx2)==localnodeidx2.AND. &
2140 & basis%NODE_POSITION_INDEX(localnodeidx1,xiidx3)==localnodeidx3)
THEN 2141 nodecount=nodecount+1
2142 nodesinline(nodecount)=localnodeidx1
2146 IF((nodecount+specialnodecount)>1)
THEN 2147 basis%NUMBER_OF_LOCAL_LINES=basis%NUMBER_OF_LOCAL_LINES+1
2148 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(basis%NUMBER_OF_LOCAL_LINES)=nodecount+specialnodecount
2149 basis%NODE_NUMBERS_IN_LOCAL_LINE(1:basis%NUMBER_OF_NODES_IN_LOCAL_LINE(basis%NUMBER_OF_LOCAL_LINES), &
2150 & basis%NUMBER_OF_LOCAL_LINES)=nodesinline(1:basis%NUMBER_OF_NODES_IN_LOCAL_LINE(basis%NUMBER_OF_LOCAL_LINES))
2151 locallineparameter=0
2152 DO locallinenodeidx=1,basis%NUMBER_OF_NODES_IN_LOCAL_LINE(basis%NUMBER_OF_LOCAL_LINES)
2153 DO derivativeidx=1,basis%NUMBER_OF_DERIVATIVES(basis%NODE_NUMBERS_IN_LOCAL_LINE(locallinenodeidx, &
2154 & basis%NUMBER_OF_LOCAL_LINES))
2155 IF(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,basis%NODE_NUMBERS_IN_LOCAL_LINE(locallinenodeidx, &
2156 & basis%NUMBER_OF_LOCAL_LINES),xiidx2)==
no_part_deriv.AND. &
2157 & basis%DERIVATIVE_ORDER_INDEX(derivativeidx,basis%NODE_NUMBERS_IN_LOCAL_LINE(locallinenodeidx, &
2159 locallineparameter=locallineparameter+1
2160 basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE(locallineparameter,basis%NUMBER_OF_LOCAL_LINES)= &
2161 & basis%ELEMENT_PARAMETER_INDEX(derivativeidx,basis%NODE_NUMBERS_IN_LOCAL_LINE(locallinenodeidx, &
2162 & basis%NUMBER_OF_LOCAL_LINES))
2163 IF(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,basis%NODE_NUMBERS_IN_LOCAL_LINE(locallinenodeidx, &
2165 basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE(locallinenodeidx,basis%NUMBER_OF_LOCAL_LINES)=derivativeidx
2170 basis%LOCAL_LINE_XI_DIRECTION(basis%NUMBER_OF_LOCAL_LINES)=xiidx1
2179 DO directionidx=-1,1,2
2186 IF(directionidx==1)
THEN 2188 localnodeidx1=maximumnodeextent(xiidx1)
2198 IF(.NOT.collapsedface)
THEN 2200 localfaceidx=localfaceidx+1
2202 DO localnodeidx3=1,maximumnodeextent(xiidx2)
2203 DO localnodeidx2=1,maximumnodeextent(xiidx3)
2205 localnodeidx=basis%NODE_POSITION_INDEX_INV(localnodeidx1,localnodeidx2,localnodeidx3,1)
2206 ELSE IF(xiidx1==2)
THEN 2207 localnodeidx=basis%NODE_POSITION_INDEX_INV(localnodeidx2,localnodeidx1,localnodeidx3,1)
2209 localnodeidx=basis%NODE_POSITION_INDEX_INV(localnodeidx2,localnodeidx3,localnodeidx1,1)
2211 IF(all(basis%NODE_NUMBERS_IN_LOCAL_FACE(1:localnodecount,localfaceidx)/=localnodeidx))
THEN 2213 localnodecount=localnodecount+1
2214 basis%NODE_NUMBERS_IN_LOCAL_FACE(localnodecount,localfaceidx)=localnodeidx
2218 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(localfaceidx)=localnodecount
2219 basis%LOCAL_FACE_XI_DIRECTION(localfaceidx)=directionidx*xiidx1
2221 localfaceparameter=0
2222 DO localnodeidx=1,basis%NUMBER_OF_NODES_IN_LOCAL_FACE(localfaceidx)
2223 localnode=basis%NODE_NUMBERS_IN_LOCAL_FACE(localnodeidx,localfaceidx)
2224 localfacederivative=0
2225 DO derivativeidx=1,basis%NUMBER_OF_DERIVATIVES(localnode)
2226 IF(basis%DERIVATIVE_ORDER_INDEX(derivativeidx,localnode,xiidx1)==
no_part_deriv)
THEN 2227 localfaceparameter=localfaceparameter+1
2228 localfacederivative=localfacederivative+1
2229 basis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE(localfacederivative,localnodeidx,localfaceidx)=derivativeidx
2230 basis%ELEMENT_PARAMETERS_IN_LOCAL_FACE(localfaceparameter,localfaceidx)= &
2231 & basis%ELEMENT_PARAMETER_INDEX(derivativeidx,localnode)
2234 basis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE(0,localnodeidx,localfaceidx)=localfacederivative
2240 CALL flagerror(
"Invalid number of xi directions.",err,error,*999)
2246 CALL flagerror(
"Basis is not a Lagrange Hermite tensor product basis.",err,error,*999)
2249 CALL flagerror(
"Basis is not associated.",err,error,*999)
2252 exits(
"Basis_LHTPBasisCreate")
2254 999
IF(
ALLOCATED(nodeatcollapse))
DEALLOCATE(nodeatcollapse)
2255 errorsexits(
"Basis_LHTPBasisCreate",err,error)
2269 INTEGER(INTG),
INTENT(IN) :: NODE_NUMBER
2270 INTEGER(INTG),
INTENT(IN) :: DERIVATIVE_NUMBER
2271 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIV_INDEX
2272 REAL(DP),
INTENT(IN) :: XI(:)
2273 INTEGER(INTG),
INTENT(OUT) :: ERR
2276 REAL(DP) :: BASIS_LHTP_BASIS_EVALUATE_DP
2278 INTEGER(INTG) :: ni,nn
2282 enters(
"BASIS_LHTP_BASIS_EVALUATE_DP",err,error,*999)
2284 basis_lhtp_basis_evaluate_dp=1.0_dp
2285 IF(
ASSOCIATED(basis))
THEN 2286 DO ni=1,basis%NUMBER_OF_XI
2287 IF(basis%NODE_AT_COLLAPSE(node_number).AND.basis%COLLAPSED_XI(ni)==
basis_xi_collapsed)
THEN 2290 SELECT CASE(basis%INTERPOLATION_TYPE(ni))
2292 SELECT CASE(basis%INTERPOLATION_ORDER(ni))
2306 local_error=
"Interpolation order value "//
trim(
number_to_vstring(basis%INTERPOLATION_ORDER(ni),
"*",err,error))// &
2308 CALL flagerror(local_error,err,error,*999)
2311 SELECT CASE(basis%INTERPOLATION_ORDER(ni))
2318 local_error=
"Interpolation order value "//
trim(
number_to_vstring(basis%INTERPOLATION_ORDER(ni),
"*",err,error))// &
2320 CALL flagerror(local_error,err,error,*999)
2324 local_error=
"Interpolation type value "//
trim(
number_to_vstring(basis%INTERPOLATION_TYPE(ni),
"*",err,error))// &
2326 CALL flagerror(local_error,err,error,*999)
2328 basis_lhtp_basis_evaluate_dp=basis_lhtp_basis_evaluate_dp*sum
2330 SELECT CASE(basis%INTERPOLATION_TYPE(ni))
2332 SELECT CASE(basis%INTERPOLATION_ORDER(ni))
2334 basis_lhtp_basis_evaluate_dp=basis_lhtp_basis_evaluate_dp* &
2338 basis_lhtp_basis_evaluate_dp=basis_lhtp_basis_evaluate_dp* &
2342 basis_lhtp_basis_evaluate_dp=basis_lhtp_basis_evaluate_dp* &
2346 local_error=
"Interpolation order value "//
trim(
number_to_vstring(basis%INTERPOLATION_ORDER(ni),
"*",err,error))// &
2348 CALL flagerror(local_error,err,error,*999)
2351 SELECT CASE(basis%INTERPOLATION_ORDER(ni))
2353 basis_lhtp_basis_evaluate_dp=basis_lhtp_basis_evaluate_dp* &
2355 & basis%DERIVATIVE_ORDER_INDEX(derivative_number,node_number,ni), &
2358 basis_lhtp_basis_evaluate_dp=basis_lhtp_basis_evaluate_dp* &
2360 & basis%DERIVATIVE_ORDER_INDEX(derivative_number,node_number,ni), &
2363 basis_lhtp_basis_evaluate_dp=basis_lhtp_basis_evaluate_dp* &
2365 & basis%DERIVATIVE_ORDER_INDEX(derivative_number,node_number,ni), &
2368 local_error=
"Interpolation order value "//
trim(
number_to_vstring(basis%INTERPOLATION_ORDER(ni),
"*",err,error))// &
2370 CALL flagerror(local_error,err,error,*999)
2374 local_error=
"Interpolation type value "//
trim(
number_to_vstring(basis%INTERPOLATION_TYPE(ni),
"*",err,error))// &
2376 CALL flagerror(local_error,err,error,*999)
2381 CALL flagerror(
"Basis is not associated",err,error,*999)
2384 exits(
"BASIS_LHTP_BASIS_EVALUATE_DP")
2386 999 errorsexits(
"BASIS_LHTP_BASIS_EVALUATE_DP",err,error)
2400 INTEGER(INTG),
INTENT(OUT) :: ERR
2403 INTEGER(INTG) :: DUMMY_ERR,ni,ni2,FACE_XI(2)
2404 LOGICAL :: LINE_BASIS_DONE,FACE_BASIS_DONE
2408 NULLIFY(new_sub_basis)
2410 enters(
"BASIS_LHTP_FAMILY_CREATE",err,error,*999)
2412 IF(
ASSOCIATED(basis))
THEN 2415 IF(basis%NUMBER_OF_XI>1)
THEN 2417 ALLOCATE(basis%LINE_BASES(basis%NUMBER_OF_XI),stat=err)
2418 IF(err/=0)
CALL flagerror(
"Could not allocate basis line bases",err,error,*999)
2419 DO ni=1,basis%NUMBER_OF_XI
2420 line_basis_done=.false.
2421 NULLIFY(new_sub_basis)
2423 IF(basis%INTERPOLATION_XI(ni2)==basis%INTERPOLATION_XI(ni).AND. &
2424 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni2)==basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni))
THEN 2425 line_basis_done=.true.
2429 IF(line_basis_done)
THEN 2430 basis%LINE_BASES(ni)%PTR=>basis%LINE_BASES(ni2)%PTR
2436 basis%LINE_BASES(ni)%PTR=>new_sub_basis
2439 IF(basis%NUMBER_OF_XI>2)
THEN 2441 ALLOCATE(basis%FACE_BASES(basis%NUMBER_OF_XI),stat=err)
2442 IF(err/=0)
CALL flagerror(
"Could not allocate basis face bases",err,error,*999)
2443 DO ni=1,basis%NUMBER_OF_XI
2447 face_basis_done=.false.
2448 NULLIFY(new_sub_basis)
2464 IF(face_basis_done)
THEN 2465 basis%FACE_BASES(ni)%PTR=>basis%FACE_BASES(ni2)%PTR
2471 new_sub_basis%LINE_BASES(1)%PTR=>basis%LINE_BASES(face_xi(1))%PTR
2472 new_sub_basis%LINE_BASES(2)%PTR=>basis%LINE_BASES(face_xi(2))%PTR
2473 basis%FACE_BASES(ni)%PTR=>new_sub_basis
2477 ALLOCATE(basis%FACE_BASES(1),stat=err)
2478 IF(err/=0)
CALL flagerror(
"Could not allocate basis face bases",err,error,*999)
2479 basis%FACE_BASES(1)%PTR=>basis
2482 ALLOCATE(basis%LINE_BASES(1),stat=err)
2483 IF(err/=0)
CALL flagerror(
"Could not allocate basis line bases",err,error,*999)
2484 basis%LINE_BASES(1)%PTR=>basis
2485 NULLIFY(basis%FACE_BASES)
2488 CALL flagerror(
"Basis is not associated",err,error,*999)
2491 exits(
"BASIS_LHTP_FAMILY_CREATE")
2493 999
IF(
ASSOCIATED(new_sub_basis))
CALL basis_family_destroy(new_sub_basis%USER_NUMBER,new_sub_basis%FAMILY_NUMBER, &
2494 & dummy_err,dummy_error,*998)
2495 998 errorsexits(
"BASIS_LHTP_FAMILY_CREATE",err,error)
2509 INTEGER(INTG),
INTENT(OUT) :: ERR
2514 enters(
"BASIS_RADIAL_FAMILY_CREATE",err,error,*999)
2516 IF(
ASSOCIATED(basis))
THEN 2518 CALL flagerror(
"Not implemented.",err,error,*999)
2520 CALL flagerror(
"Basis is not associated",err,error,*999)
2523 exits(
"BASIS_RADIAL_FAMILY_CREATE")
2525 999 errorsexits(
"BASIS_RADIAL_FAMILY_CREATE",err,error)
2539 INTEGER(INTG),
INTENT(IN) :: LOCAL_NODE_NUMBER
2540 REAL(DP),
INTENT(OUT) :: XI(:)
2541 INTEGER(INTG),
INTENT(OUT) :: ERR
2544 INTEGER(INTG) :: xi_idx
2547 enters(
"BASIS_LOCAL_NODE_XI_CALCULATE",err,error,*999)
2549 IF(
ASSOCIATED(basis))
THEN 2550 IF(basis%BASIS_FINISHED)
THEN 2551 IF(local_node_number>0.AND.local_node_number<=basis%NUMBER_OF_NODES)
THEN 2552 IF(
SIZE(xi,1)>=basis%NUMBER_OF_XI)
THEN 2553 SELECT CASE(basis%TYPE)
2555 DO xi_idx=1,basis%NUMBER_OF_XI
2556 xi(xi_idx)=
REAL(basis%node_position_index(local_node_number,xi_idx)-1,
dp)/ &
2557 & REAL(BASIS%NUMBER_OF_NODES_XIC(xi_idx)-1,DP)
2560 DO xi_idx=1,basis%NUMBER_OF_XI
2561 xi(xi_idx)=
REAL(BASIS%NUMBER_OF_NODES_XIC(xi_idx)-BASIS%NODE_POSITION_INDEX(LOCAL_NODE_NUMBER,xi_idx),DP)/ &
2562 & REAL(BASIS%NUMBER_OF_NODES_XIC(xi_idx)-1,DP)
2565 CALL flagerror(
"Not implemented.",err,error,*999)
2567 CALL flagerror(
"Not implemented.",err,error,*999)
2569 CALL flagerror(
"Not implemented.",err,error,*999)
2571 CALL flagerror(
"Not implemented.",err,error,*999)
2573 CALL flagerror(
"Not implemented.",err,error,*999)
2577 CALL flagerror(local_error,err,error,*999)
2580 local_error=
"The size of the specified xic array of "//
trim(
number_to_vstring(
SIZE(xi,1),
"*",err,error))// &
2581 &
" is invalid. The size of the xi array must be >= "// &
2583 CALL flagerror(local_error,err,error,*999)
2586 local_error=
"The specified local node number of "//
trim(
number_to_vstring(local_node_number,
"*",err,error))// &
2587 &
" is invalid. The local node number must be > 0 and <= "// &
2589 CALL flagerror(local_error,err,error,*999)
2592 CALL flagerror(
"Basis has not been finished.",err,error,*999)
2595 CALL flagerror(
"Basis is not associated",err,error,*999)
2598 exits(
"BASIS_LOCAL_NODE_XI_CALCULATE")
2600 999 errorsexits(
"BASIS_LOCAL_NODE_XI_CALCULATE",err,error)
2613 INTEGER(INTG),
INTENT(OUT) :: NUMBER_OF_LOCAL_NODES
2614 INTEGER(INTG),
INTENT(OUT) :: ERR
2618 enters(
"BASIS_NUMBER_OF_LOCAL_NODES_GET",err,error,*999)
2620 IF(
ASSOCIATED(basis))
THEN 2621 number_of_local_nodes=basis%NUMBER_OF_NODES
2623 CALL flagerror(
"Basis is not associated",err,error,*999)
2626 exits(
"BASIS_NUMBER_OF_LOCAL_NODES_GET")
2628 999 errorsexits(
"BASIS_NUMBER_OF_LOCAL_NODES_GET",err,error)
2641 INTEGER(INTG),
INTENT(OUT) :: NUMBER_OF_XI
2642 INTEGER(INTG),
INTENT(OUT) :: ERR
2646 enters(
"BASIS_NUMBER_OF_XI_GET",err,error,*999)
2648 IF(
ASSOCIATED(basis))
THEN 2649 IF(basis%BASIS_FINISHED)
THEN 2650 number_of_xi=basis%NUMBER_OF_XI
2652 CALL flagerror(
"Basis has not been finished.",err,error,*999)
2655 CALL flagerror(
"Basis is not associated.",err,error,*999)
2658 exits(
"BASIS_NUMBER_OF_XI_GET")
2660 999 errorsexits(
"BASIS_NUMBER_OF_XI_GET",err,error)
2672 INTEGER(INTG),
INTENT(IN) :: USER_NUMBER
2673 INTEGER(INTG),
INTENT(IN) :: NUMBER_OF_XI
2674 INTEGER(INTG),
INTENT(OUT) :: ERR
2679 enters(
"BASIS_NUMBER_OF_XI_SET_NUMBER",err,error,*999)
2684 exits(
"BASIS_NUMBER_OF_XI_SET_NUMBER")
2686 999 errorsexits(
"BASIS_NUMBER_OF_XI_SET_NUMBER",err,error)
2699 INTEGER(INTG),
INTENT(IN) :: NUMBER_OF_XI
2700 INTEGER(INTG),
INTENT(OUT) :: ERR
2703 INTEGER(INTG) :: OLD_INTERPOLATION_XI(3),OLD_NUMBER_OF_GAUSS_XI(3),OLD_COLLAPSED_XI(3)
2706 enters(
"BASIS_NUMBER_OF_XI_SET_PTR",err,error,*999)
2708 IF(
ASSOCIATED(basis))
THEN 2709 IF(basis%BASIS_FINISHED)
THEN 2710 CALL flagerror(
"Basis has been finished",err,error,*999)
2712 SELECT CASE(basis%TYPE)
2714 IF(number_of_xi>0.AND.number_of_xi<4)
THEN 2715 IF(basis%NUMBER_OF_XI/=number_of_xi)
THEN 2717 old_interpolation_xi=basis%INTERPOLATION_XI
2718 old_collapsed_xi=basis%COLLAPSED_XI
2719 DEALLOCATE(basis%INTERPOLATION_XI)
2720 DEALLOCATE(basis%COLLAPSED_XI)
2721 ALLOCATE(basis%INTERPOLATION_XI(number_of_xi),stat=err)
2722 IF(err/=0)
CALL flagerror(
"Could not allocate interpolation type",err,error,*999)
2723 ALLOCATE(basis%COLLAPSED_XI(number_of_xi),stat=err)
2724 IF(err/=0)
CALL flagerror(
"Could not allocate collapsed xi",err,error,*999)
2725 IF(number_of_xi>basis%NUMBER_OF_XI)
THEN 2726 basis%INTERPOLATION_XI(1:basis%NUMBER_OF_XI)=old_interpolation_xi(1:basis%NUMBER_OF_XI)
2727 basis%INTERPOLATION_XI(basis%NUMBER_OF_XI+1:number_of_xi)=old_interpolation_xi(1)
2728 basis%COLLAPSED_XI(1:basis%NUMBER_OF_XI)=old_collapsed_xi(1:basis%NUMBER_OF_XI)
2729 basis%COLLAPSED_XI(basis%NUMBER_OF_XI+1:number_of_xi)=old_collapsed_xi(1)
2731 basis%INTERPOLATION_XI(1:number_of_xi)=old_interpolation_xi(1:number_of_xi)
2732 basis%COLLAPSED_XI(1:number_of_xi)=old_collapsed_xi(1:number_of_xi)
2735 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 2736 old_number_of_gauss_xi=basis%QUADRATURE%NUMBER_OF_GAUSS_XI
2737 DEALLOCATE(basis%QUADRATURE%NUMBER_OF_GAUSS_XI)
2738 ALLOCATE(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(number_of_xi),stat=err)
2739 IF(err/=0)
CALL flagerror(
"Could not allocate number of Gauss xi",err,error,*999)
2740 IF(number_of_xi>basis%NUMBER_OF_XI)
THEN 2741 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1:basis%NUMBER_OF_XI)=old_number_of_gauss_xi(1:basis%NUMBER_OF_XI)
2742 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(basis%NUMBER_OF_XI+1:number_of_xi)=old_number_of_gauss_xi(1)
2744 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1:number_of_xi)=old_number_of_gauss_xi(1:number_of_xi)
2747 basis%NUMBER_OF_XI=number_of_xi
2750 local_error=
"Invalid number of xi directions specified ("// &
2752 &
") for a Lagrange-Hermite basis. You must specify between 1 and 3 xi directions" 2753 CALL flagerror(local_error,err,error,*999)
2756 IF(number_of_xi>1.AND.number_of_xi<4)
THEN 2757 IF(basis%NUMBER_OF_XI/=number_of_xi)
THEN 2759 old_interpolation_xi=basis%INTERPOLATION_XI
2760 old_collapsed_xi=basis%COLLAPSED_XI
2761 DEALLOCATE(basis%INTERPOLATION_XI)
2762 DEALLOCATE(basis%COLLAPSED_XI)
2763 ALLOCATE(basis%INTERPOLATION_XI(number_of_xi),stat=err)
2764 IF(err/=0)
CALL flagerror(
"Could not allocate interpolation type",err,error,*999)
2765 ALLOCATE(basis%COLLAPSED_XI(number_of_xi),stat=err)
2766 IF(err/=0)
CALL flagerror(
"Could not allocate collapsed xi.",err,error,*999)
2767 IF(number_of_xi>basis%NUMBER_OF_XI)
THEN 2768 basis%INTERPOLATION_XI(1:basis%NUMBER_OF_XI)=old_interpolation_xi(1:basis%NUMBER_OF_XI)
2769 basis%INTERPOLATION_XI(basis%NUMBER_OF_XI+1:number_of_xi)=old_interpolation_xi(1)
2770 basis%COLLAPSED_XI(1:basis%NUMBER_OF_XI)=old_collapsed_xi(1:basis%NUMBER_OF_XI)
2771 basis%COLLAPSED_XI(basis%NUMBER_OF_XI+1:number_of_xi)=old_collapsed_xi(1)
2773 basis%INTERPOLATION_XI(1:number_of_xi)=old_interpolation_xi(1:number_of_xi)
2774 basis%COLLAPSED_XI(1:number_of_xi)=old_collapsed_xi(1:number_of_xi)
2776 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 2777 old_number_of_gauss_xi=basis%QUADRATURE%NUMBER_OF_GAUSS_XI
2778 DEALLOCATE(basis%QUADRATURE%NUMBER_OF_GAUSS_XI)
2779 ALLOCATE(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(number_of_xi),stat=err)
2780 IF(err/=0)
CALL flagerror(
"Could not allocate number of Gauss xi",err,error,*999)
2781 IF(number_of_xi>basis%NUMBER_OF_XI)
THEN 2782 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1:basis%NUMBER_OF_XI)=old_number_of_gauss_xi(1:basis%NUMBER_OF_XI)
2783 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(basis%NUMBER_OF_XI+1:number_of_xi)=old_number_of_gauss_xi(1)
2785 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1:number_of_xi)=old_number_of_gauss_xi(1:number_of_xi)
2788 basis%NUMBER_OF_XI=number_of_xi
2791 local_error=
"Invalid number of xi directions specified ("// &
2793 &
") for a simplex basis. You must specify between 2 and 3 xi directions" 2794 CALL flagerror(local_error,err,error,*999)
2797 CALL flagerror(
"Basis type invalid or not implemented",err,error,*999)
2801 CALL flagerror(
"Basis is not associated",err,error,*999)
2804 exits(
"BASIS_NUMBER_OF_XI_SET_PTR")
2806 999 errorsexits(
"BASIS_NUMBER_OF_XI_SET_PTR",err,error)
2819 INTEGER(INTG),
INTENT(OUT) :: ERR
2822 INTEGER(INTG) :: scheme_idx,i,j,k,MAX_NUM_GAUSS,ng,ni,nk,nn,ns,nu,NUM_GAUSS_1,NUM_GAUSS_2,NUM_GAUSS_3
2823 REAL(DP) :: XI(3),GSX(4,20),GSW(20)
2824 REAL(DP),
ALLOCATABLE :: POSITIONS(:,:),POSITIONS_MATRIX(:,:,:,:),WEIGHTS(:,:)
2828 INTEGER(INTG) :: MAX_NUM_FACE_GAUSS,face_idx,NORMAL,FACE_XI(2),numberOfFaceXiCoordinates
2831 NULLIFY(new_schemes)
2833 enters(
"BASIS_QUADRATURE_CREATE",err,error,*999)
2835 IF(
ASSOCIATED(basis))
THEN 2836 IF(
ASSOCIATED(basis%QUADRATURE%SCHEMES))
THEN 2837 local_error=
"The quadrature schemes on basis number "//
trim(
number_to_vstring(basis%USER_NUMBER,
"*",err,error))// &
2838 &
" are already associated" 2839 CALL flagerror(local_error,err,error,*998)
2846 SELECT CASE(basis%QUADRATURE%TYPE)
2849 ALLOCATE(new_scheme,stat=err)
2850 IF(err/=0)
CALL flagerror(
"Could not allocate new quadrature scheme",err,error,*999)
2851 new_scheme%QUADRATURE=>basis%QUADRATURE
2852 basis%QUADRATURE%NUMBER_OF_SCHEMES=1
2853 ALLOCATE(new_schemes(basis%QUADRATURE%NUMBER_OF_SCHEMES),stat=err)
2854 IF(err/=0)
CALL flagerror(
"Could not allocate new quadratures scheme",err,error,*999)
2855 new_schemes(1)%PTR=>new_scheme
2856 basis%QUADRATURE%SCHEMES=>new_schemes
2859 IF(err/=0)
CALL flagerror(
"Could not allocate quadrature scheme map",err,error,*999)
2861 NULLIFY(basis%QUADRATURE%QUADRATURE_SCHEME_MAP(scheme_idx)%PTR)
2865 new_scheme%NUMBER_OF_GAUSS=1
2867 DO ni=1,basis%NUMBER_OF_XI
2868 new_scheme%NUMBER_OF_GAUSS=new_scheme%NUMBER_OF_GAUSS*basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)
2869 IF(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)>max_num_gauss) max_num_gauss=basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)
2871 ALLOCATE(new_scheme%GAUSS_POSITIONS(basis%NUMBER_OF_XI_COORDINATES,new_scheme%NUMBER_OF_GAUSS),stat=err)
2872 IF(err/=0)
CALL flagerror(
"Could not allocate Gauss positions",err,error,*999)
2873 ALLOCATE(new_scheme%GAUSS_WEIGHTS(new_scheme%NUMBER_OF_GAUSS),stat=err)
2874 IF(err/=0)
CALL flagerror(
"Could not allocate Gauss weights",err,error,*999)
2875 ALLOCATE(new_scheme%GAUSS_BASIS_FNS(basis%NUMBER_OF_ELEMENT_PARAMETERS,basis%NUMBER_OF_PARTIAL_DERIVATIVES, &
2876 & new_scheme%NUMBER_OF_GAUSS),stat=err)
2877 IF(err/=0)
CALL flagerror(
"Could not allocate Gauss basis functions",err,error,*999)
2878 ALLOCATE(weights(max_num_gauss,3),stat=err)
2879 IF(err/=0)
CALL flagerror(
"Could not allocate weights",err,error,*999)
2880 ALLOCATE(positions(max_num_gauss,3),stat=err)
2881 IF(err/=0)
CALL flagerror(
"Could not allocate positions",err,error,*999)
2882 ALLOCATE(positions_matrix(max_num_gauss,max_num_gauss,max_num_gauss,3),stat=err)
2883 IF(err/=0)
CALL flagerror(
"Could not allocate positions matrix",err,error,*999)
2886 positions_matrix=0.0_dp
2887 DO ni=1,basis%NUMBER_OF_XI
2888 CALL gauss_legendre(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni),0.0_dp,1.0_dp, &
2889 & positions(1:basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni),ni), &
2890 & weights(1:basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni),ni),err,error,*999)
2892 SELECT CASE(basis%NUMBER_OF_XI)
2894 num_gauss_1=basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1)
2898 num_gauss_1=basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1)
2899 num_gauss_2=basis%QUADRATURE%NUMBER_OF_GAUSS_XI(2)
2902 num_gauss_1=basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1)
2903 num_gauss_2=basis%QUADRATURE%NUMBER_OF_GAUSS_XI(2)
2904 num_gauss_3=basis%QUADRATURE%NUMBER_OF_GAUSS_XI(3)
2906 CALL flagerror(
"Invalid number of xi directions",err,error,*999)
2911 positions_matrix(i,j,k,1)=positions(i,1)
2912 positions_matrix(i,j,k,2)=positions(j,2)
2913 positions_matrix(i,j,k,3)=positions(k,3)
2914 xi(1:basis%NUMBER_OF_XI)=positions_matrix(i,j,k,1:basis%NUMBER_OF_XI)
2915 ng=i+(j-1+(k-1)*num_gauss_2)*num_gauss_1
2916 new_scheme%GAUSS_WEIGHTS(ng)=weights(i,1)*weights(j,2)*weights(k,3)
2917 new_scheme%GAUSS_POSITIONS(1:basis%NUMBER_OF_XI_COORDINATES,ng)=xi(1:basis%NUMBER_OF_XI_COORDINATES)
2919 DO nn=1,basis%NUMBER_OF_NODES
2920 DO nk=1,basis%NUMBER_OF_DERIVATIVES(nn)
2922 DO nu=1,basis%NUMBER_OF_PARTIAL_DERIVATIVES
2923 SELECT CASE(basis%TYPE)
2928 CALL flagerror(
"Not implemented",err,error,*999)
2937 IF(basis%QUADRATURE%EVALUATE_FACE_GAUSS)
THEN 2938 IF(basis%NUMBER_OF_XI==3)
THEN 2940 max_num_face_gauss=product(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1:basis%NUMBER_OF_XI))
2941 max_num_face_gauss=max_num_face_gauss/minval(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1:basis%NUMBER_OF_XI))
2942 ALLOCATE(new_scheme%NUMBER_OF_FACE_GAUSS(basis%NUMBER_OF_LOCAL_FACES),stat=err)
2943 IF(err/=0)
CALL flagerror(
"Could not allocate number of face gauss",err,error,*999)
2944 ALLOCATE(new_scheme%FACE_GAUSS_POSITIONS(basis%NUMBER_OF_XI_COORDINATES,max_num_face_gauss, &
2945 & basis%NUMBER_OF_LOCAL_FACES),stat=err)
2946 IF(err/=0)
CALL flagerror(
"Could not allocate face Gauss positions",err,error,*999)
2947 ALLOCATE(new_scheme%FACE_GAUSS_WEIGHTS(max_num_face_gauss,basis%NUMBER_OF_LOCAL_FACES),stat=err)
2948 IF(err/=0)
CALL flagerror(
"Could not allocate face Gauss weights",err,error,*999)
2949 ALLOCATE(new_scheme%FACE_GAUSS_BASIS_FNS(basis%NUMBER_OF_ELEMENT_PARAMETERS,basis%NUMBER_OF_PARTIAL_DERIVATIVES, &
2950 & max_num_face_gauss,basis%NUMBER_OF_LOCAL_FACES),stat=err)
2951 IF(err/=0)
CALL flagerror(
"Could not allocate face Gauss basis function values array",err,error,*999)
2953 new_scheme%FACE_GAUSS_POSITIONS=0.0_dp
2954 new_scheme%FACE_GAUSS_WEIGHTS=0.0_dp
2955 new_scheme%FACE_GAUSS_BASIS_FNS=0.0_dp
2957 DO face_idx=1,basis%NUMBER_OF_LOCAL_FACES
2959 normal=basis%LOCAL_FACE_XI_DIRECTION(face_idx)
2960 IF(normal<0_intg)
THEN 2961 xi(abs(normal))=0.0_dp
2963 xi(abs(normal))=1.0_dp
2968 new_scheme%NUMBER_OF_FACE_GAUSS(face_idx)=product(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(face_xi))
2970 DO j=1,basis%QUADRATURE%NUMBER_OF_GAUSS_XI(face_xi(2))
2971 xi(face_xi(2))=positions(j,face_xi(2))
2972 DO i=1,basis%QUADRATURE%NUMBER_OF_GAUSS_XI(face_xi(1))
2973 xi(face_xi(1))=positions(i,face_xi(1))
2976 new_scheme%FACE_GAUSS_WEIGHTS(ng,face_idx)=weights(i,face_xi(1))*weights(j,face_xi(2))
2977 new_scheme%FACE_GAUSS_POSITIONS(1:3,ng,face_idx)=xi(1:3)
2980 DO nn=1,basis%NUMBER_OF_NODES
2981 DO nk=1,basis%NUMBER_OF_DERIVATIVES(nn)
2983 DO nu=1,basis%NUMBER_OF_PARTIAL_DERIVATIVES
2984 SELECT CASE(basis%TYPE)
2986 new_scheme%FACE_GAUSS_BASIS_FNS(ns,nu,ng,face_idx)= &
2990 CALL flagerror(
"Not implemented",err,error,*999)
3000 CALL flagerror(
"Cannot evaluate face quadrature schemes for a non three dimensional element.",err,error,*999)
3005 DEALLOCATE(positions)
3006 DEALLOCATE(positions_matrix)
3008 CALL flagerror(
"Gauss Laguerre quadrature type not implemented.",err,error,*999)
3010 CALL flagerror(
"Gauss Hermite quadrature type not implemented.",err,error,*999)
3013 ALLOCATE(new_scheme,stat=err)
3014 IF(err/=0)
CALL flagerror(
"Could not allocate new quadrature scheme.",err,error,*999)
3015 new_scheme%QUADRATURE=>basis%QUADRATURE
3016 basis%QUADRATURE%NUMBER_OF_SCHEMES=1
3017 ALLOCATE(new_schemes(basis%QUADRATURE%NUMBER_OF_SCHEMES),stat=err)
3018 IF(err/=0)
CALL flagerror(
"Could not allocate new quadratures scheme.",err,error,*999)
3019 new_schemes(1)%PTR=>new_scheme
3020 basis%QUADRATURE%SCHEMES=>new_schemes
3023 IF(err/=0)
CALL flagerror(
"Could not allocate quadrature scheme map.",err,error,*999)
3025 NULLIFY(basis%QUADRATURE%QUADRATURE_SCHEME_MAP(scheme_idx)%PTR)
3029 CALL gauss_simplex(basis%QUADRATURE%GAUSS_ORDER,basis%NUMBER_OF_XI_COORDINATES,new_scheme%NUMBER_OF_GAUSS,gsx,gsw, &
3031 ALLOCATE(new_scheme%GAUSS_POSITIONS(basis%NUMBER_OF_XI_COORDINATES,new_scheme%NUMBER_OF_GAUSS),stat=err)
3032 IF(err/=0)
CALL flagerror(
"Could not allocate Gauss positions.",err,error,*999)
3033 ALLOCATE(new_scheme%GAUSS_WEIGHTS(new_scheme%NUMBER_OF_GAUSS),stat=err)
3034 IF(err/=0)
CALL flagerror(
"Could not allocate Gauss weights.",err,error,*999)
3035 ALLOCATE(new_scheme%GAUSS_BASIS_FNS(basis%NUMBER_OF_ELEMENT_PARAMETERS,basis%NUMBER_OF_PARTIAL_DERIVATIVES, &
3036 & new_scheme%NUMBER_OF_GAUSS),stat=err)
3037 IF(err/=0)
CALL flagerror(
"Could not allocate Gauss basis functions.",err,error,*999)
3038 new_scheme%GAUSS_POSITIONS(1:basis%NUMBER_OF_XI_COORDINATES,1:new_scheme%NUMBER_OF_GAUSS)= &
3039 & gsx(1:basis%NUMBER_OF_XI_COORDINATES,1:new_scheme%NUMBER_OF_GAUSS)
3040 new_scheme%GAUSS_WEIGHTS(1:new_scheme%NUMBER_OF_GAUSS)=gsw(1:new_scheme%NUMBER_OF_GAUSS)
3041 DO ng=1,new_scheme%NUMBER_OF_GAUSS
3043 DO nn=1,basis%NUMBER_OF_NODES
3044 DO nk=1,basis%NUMBER_OF_DERIVATIVES(nn)
3046 DO nu=1,basis%NUMBER_OF_PARTIAL_DERIVATIVES
3047 SELECT CASE(basis%TYPE)
3050 new_scheme%GAUSS_BASIS_FNS(ns,nu,ng)= &
3055 CALL flagerror(
"Not implemented.",err,error,*999)
3062 IF(basis%QUADRATURE%EVALUATE_FACE_GAUSS)
THEN 3063 IF(basis%NUMBER_OF_XI==3)
THEN 3065 max_num_face_gauss=product(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1:basis%NUMBER_OF_XI))
3066 max_num_face_gauss=max_num_face_gauss/minval(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1:basis%NUMBER_OF_XI))
3067 ALLOCATE(new_scheme%NUMBER_OF_FACE_GAUSS(basis%NUMBER_OF_LOCAL_FACES),stat=err)
3068 IF(err/=0)
CALL flagerror(
"Could not allocate number of face gauss",err,error,*999)
3069 ALLOCATE(new_scheme%FACE_GAUSS_POSITIONS(basis%NUMBER_OF_XI_COORDINATES,max_num_face_gauss, &
3070 & basis%NUMBER_OF_LOCAL_FACES),stat=err)
3071 IF(err/=0)
CALL flagerror(
"Could not allocate face Gauss positions",err,error,*999)
3072 ALLOCATE(new_scheme%FACE_GAUSS_WEIGHTS(max_num_face_gauss,basis%NUMBER_OF_LOCAL_FACES),stat=err)
3073 IF(err/=0)
CALL flagerror(
"Could not allocate face Gauss weights",err,error,*999)
3074 ALLOCATE(new_scheme%FACE_GAUSS_BASIS_FNS(basis%NUMBER_OF_ELEMENT_PARAMETERS,basis%NUMBER_OF_PARTIAL_DERIVATIVES, &
3075 & max_num_face_gauss,basis%NUMBER_OF_LOCAL_FACES),stat=err)
3076 IF(err/=0)
CALL flagerror(
"Could not allocate face Gauss basis function values array",err,error,*999)
3078 new_scheme%FACE_GAUSS_POSITIONS=0.0_dp
3079 new_scheme%FACE_GAUSS_WEIGHTS=0.0_dp
3080 new_scheme%FACE_GAUSS_BASIS_FNS=0.0_dp
3082 DO face_idx=1,basis%NUMBER_OF_LOCAL_FACES
3084 numberoffacexicoordinates = basis%NUMBER_OF_XI
3086 CALL gauss_simplex(basis%QUADRATURE%GAUSS_ORDER,numberoffacexicoordinates, &
3087 & new_scheme%NUMBER_OF_FACE_GAUSS(face_idx),gsx,gsw,err,error,*999)
3088 IF(err/=0)
CALL flagerror(
"Could not allocate Gauss basis functions",err,error,*999)
3089 new_scheme%FACE_GAUSS_POSITIONS(1:numberoffacexicoordinates,1:new_scheme%NUMBER_OF_FACE_GAUSS(face_idx), &
3090 & face_idx)=gsx(1:numberoffacexicoordinates,1:new_scheme%NUMBER_OF_FACE_GAUSS(face_idx))
3091 new_scheme%FACE_GAUSS_WEIGHTS(1:new_scheme%NUMBER_OF_FACE_GAUSS(face_idx),face_idx)= &
3092 & gsw(1:new_scheme%NUMBER_OF_FACE_GAUSS(face_idx))
3094 DO ng=1,new_scheme%NUMBER_OF_FACE_GAUSS(face_idx)
3096 DO nn=1,basis%NUMBER_OF_NODES_IN_LOCAL_FACE(face_idx)
3097 DO nk=1,basis%NUMBER_OF_DERIVATIVES(nn)
3099 DO nu=1,basis%NUMBER_OF_PARTIAL_DERIVATIVES
3100 SELECT CASE(basis%TYPE)
3102 new_scheme%FACE_GAUSS_BASIS_FNS(ns,nu,ng,face_idx)= &
3104 & new_scheme%FACE_GAUSS_POSITIONS(1:numberoffacexicoordinates,ng,face_idx),err,error)
3107 CALL flagerror(
"Not implemented",err,error,*999)
3116 CALL flagerror(
"Cannot evaluate face quadrature schemes for a non three dimensional element.",err,error,*999)
3120 local_error=
"Quadrature type "//
trim(
number_to_vstring(basis%QUADRATURE%TYPE,
"*",err,error))//
" is invalid." 3121 CALL flagerror(local_error,err,error,*999)
3125 CALL flagerror(
"Basis is not associated.",err,error,*998)
3131 &
'(" Number of gauss points(ni):",3(X,I2))',
'(22X,3(X,I2))',err,error,*999)
3134 DO scheme_idx=1,basis%QUADRATURE%NUMBER_OF_SCHEMES
3135 scheme=>basis%QUADRATURE%SCHEMES(scheme_idx)%PTR
3141 DO ng=1,scheme%NUMBER_OF_GAUSS
3144 &
'(" position(ni) :",3(X,F12.4))',
'(26X,3(X,F12.4))',err,error,*999)
3146 &
"(F12.4)",err,error,*999)
3151 DO ng=1,scheme%NUMBER_OF_GAUSS
3153 DO nu=1,basis%NUMBER_OF_PARTIAL_DERIVATIVES
3156 & scheme%GAUSS_BASIS_FNS(:,nu,ng),
'(" BASIS FNS(ns) :",4(X,F12.4))',
'(26X,4(X,F12.4))',err,error,*999)
3163 exits(
"BASIS_QUADRATURE_CREATE")
3165 999
IF(
ASSOCIATED(new_scheme))
THEN 3166 IF(
ALLOCATED(new_scheme%GAUSS_POSITIONS))
DEALLOCATE(new_scheme%GAUSS_POSITIONS)
3167 IF(
ALLOCATED(new_scheme%GAUSS_WEIGHTS))
DEALLOCATE(new_scheme%GAUSS_WEIGHTS)
3168 IF(
ALLOCATED(new_scheme%GAUSS_BASIS_FNS))
DEALLOCATE(new_scheme%GAUSS_BASIS_FNS)
3169 DEALLOCATE(new_scheme)
3171 IF(
ALLOCATED(weights))
DEALLOCATE(weights)
3172 IF(
ALLOCATED(positions))
DEALLOCATE(positions)
3173 IF(
ALLOCATED(positions_matrix))
DEALLOCATE(positions_matrix)
3174 IF(
ASSOCIATED(new_schemes))
DEALLOCATE(new_schemes)
3175 NULLIFY(basis%QUADRATURE%SCHEMES)
3176 998 errorsexits(
"BASIS_QUADRATURE_CREATE",err,error)
3190 INTEGER(INTG),
INTENT(OUT) :: ERR
3194 enters(
"BASIS_QUADRATURE_DESTROY",err,error,*999)
3196 IF(
ASSOCIATED(quadrature))
THEN 3197 CALL flagerror(
"Not implemented.",err,error,*999)
3199 CALL flagerror(
"Basis is not associated",err,error,*999)
3202 exits(
"BASIS_QUADRATURE_DESTROY")
3204 999 errorsexits(
"BASIS_QUADRATURE_DESTROY",err,error)
3217 INTEGER(INTG),
INTENT(OUT) :: ERR
3220 INTEGER(INTG) :: scheme_idx
3223 enters(
"BASIS_QUADRATURE_FINALISE",err,error,*999)
3225 IF(
ASSOCIATED(basis))
THEN 3226 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 3227 DO scheme_idx=1,basis%QUADRATURE%NUMBER_OF_SCHEMES
3228 scheme=>basis%QUADRATURE%SCHEMES(scheme_idx)%PTR
3230 IF (
ASSOCIATED(scheme))
THEN 3231 IF(
ALLOCATED(scheme%GAUSS_POSITIONS))
DEALLOCATE(scheme%GAUSS_POSITIONS)
3232 IF(
ALLOCATED(scheme%GAUSS_WEIGHTS))
DEALLOCATE(scheme%GAUSS_WEIGHTS)
3233 IF(
ALLOCATED(scheme%GAUSS_BASIS_FNS))
DEALLOCATE(scheme%GAUSS_BASIS_FNS)
3237 IF(
ASSOCIATED(basis%QUADRATURE%SCHEMES))
DEALLOCATE(basis%QUADRATURE%SCHEMES)
3238 basis%QUADRATURE%NUMBER_OF_SCHEMES=0
3239 IF(
ALLOCATED(basis%QUADRATURE%NUMBER_OF_GAUSS_XI))
DEALLOCATE(basis%QUADRATURE%NUMBER_OF_GAUSS_XI)
3240 NULLIFY(basis%QUADRATURE%BASIS)
3241 basis%QUADRATURE%TYPE=-1
3242 IF(
ALLOCATED(basis%QUADRATURE%QUADRATURE_SCHEME_MAP))
DEALLOCATE(basis%QUADRATURE%QUADRATURE_SCHEME_MAP)
3244 CALL flagerror(
"Basis quadrature basis is not associated",err,error,*999)
3247 CALL flagerror(
"Basis is not associated",err,error,*999)
3250 exits(
"BASIS_QUADRATURE_FINALISE")
3252 999 errorsexits(
"BASIS_QUADRATURE_FINALISE",err,error)
3265 INTEGER(INTG),
INTENT(OUT) :: ERR
3271 enters(
"BASIS_QUADRATURE_INITIALISE",err,error,*999)
3273 IF(
ASSOCIATED(basis))
THEN 3274 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 3276 &
" already has a quadrature associated" 3277 CALL flagerror(local_error,err,error,*998)
3279 SELECT CASE(basis%TYPE)
3282 basis%QUADRATURE%NUMBER_OF_SCHEMES=0
3283 NULLIFY(basis%QUADRATURE%SCHEMES)
3284 basis%QUADRATURE%BASIS=>basis
3287 ALLOCATE(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(basis%NUMBER_OF_XI),stat=err)
3288 IF(err/=0)
CALL flagerror(
"Could not allocate number of Gauss in each xi direction",err,error,*999)
3289 DO ni=1,basis%NUMBER_OF_XI
3290 SELECT CASE(basis%INTERPOLATION_XI(ni))
3292 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)=2
3294 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)=3
3296 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)=4
3298 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)=3
3300 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)=4
3302 local_error=
"Interpolation xi value "//
trim(
number_to_vstring(basis%INTERPOLATION_XI(ni),
"*",err,error))// &
3304 CALL flagerror(local_error,err,error,*999)
3309 basis%QUADRATURE%NUMBER_OF_SCHEMES=0
3310 NULLIFY(basis%QUADRATURE%SCHEMES)
3311 basis%QUADRATURE%BASIS=>basis
3314 ALLOCATE(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(basis%NUMBER_OF_XI),stat=err)
3315 IF(err/=0)
CALL flagerror(
"Could not allocate number of Gauss in each xi direction",err,error,*999)
3317 SELECT CASE(basis%INTERPOLATION_XI(1))
3319 SELECT CASE(basis%NUMBER_OF_XI)
3321 basis%QUADRATURE%GAUSS_ORDER=2
3323 basis%QUADRATURE%GAUSS_ORDER=3
3325 basis%QUADRATURE%GAUSS_ORDER=3
3327 local_error=
"The number of xi directions ("//
trim(
number_to_vstring(basis%NUMBER_OF_XI,
"*",err,error))// &
3329 CALL flagerror(local_error,err,error,*999)
3332 SELECT CASE(basis%NUMBER_OF_XI)
3334 basis%QUADRATURE%GAUSS_ORDER=3
3336 basis%QUADRATURE%GAUSS_ORDER=3
3338 basis%QUADRATURE%GAUSS_ORDER=5
3340 local_error=
"The number of xi directions ("//
trim(
number_to_vstring(basis%NUMBER_OF_XI,
"*",err,error))// &
3342 CALL flagerror(local_error,err,error,*999)
3345 SELECT CASE(basis%NUMBER_OF_XI)
3347 basis%QUADRATURE%GAUSS_ORDER=3
3349 basis%QUADRATURE%GAUSS_ORDER=5
3351 basis%QUADRATURE%GAUSS_ORDER=5
3353 local_error=
"The number of xi directions ("//
trim(
number_to_vstring(basis%NUMBER_OF_XI,
"*",err,error))// &
3355 CALL flagerror(local_error,err,error,*999)
3358 local_error=
"Interpolation xi value "//
trim(
number_to_vstring(basis%INTERPOLATION_XI(1),
"*",err,error))// &
3359 &
" in xi direction 1 is invalid" 3360 CALL flagerror(local_error,err,error,*999)
3362 basis%QUADRATURE%NUMBER_OF_GAUSS_XI=basis%QUADRATURE%GAUSS_ORDER
3364 local_error=
"Basis type value "//
trim(
number_to_vstring(basis%INTERPOLATION_XI(ni),
"*",err,error))// &
3365 &
" is invalid or not implemented" 3366 CALL flagerror(local_error,err,error,*999)
3370 CALL flagerror(
"Basis is not associated",err,error,*998)
3373 exits(
"BASIS_QUADRATURE_INITIALISE")
3375 999
IF(
ALLOCATED(basis%QUADRATURE%NUMBER_OF_GAUSS_XI))
DEALLOCATE(basis%QUADRATURE%NUMBER_OF_GAUSS_XI)
3376 998 errorsexits(
"BASIS_QUADRATURE_INITIALISE",err,error)
3390 INTEGER(INTG),
INTENT(OUT) :: QUADRATURE_NUMBER_OF_GAUSS_XI(:)
3391 INTEGER(INTG),
INTENT(OUT) :: ERR
3396 enters(
"BASIS_QUADRATURE_NUMBER_OF_GAUSS_XI_GET",err,error,*999)
3398 IF(
ASSOCIATED(basis))
THEN 3399 IF(basis%BASIS_FINISHED)
THEN 3400 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 3401 IF(
SIZE(quadrature_number_of_gauss_xi,1)>=
SIZE(basis%QUADRATURE%NUMBER_OF_GAUSS_XI,1))
THEN 3402 quadrature_number_of_gauss_xi=basis%QUADRATURE%NUMBER_OF_GAUSS_XI
3404 local_error=
"The size of QUADRATURE_NUMBER_OF_GAUSS_XI is too small. The supplied size is "// &
3405 &
trim(
number_to_vstring(
SIZE(quadrature_number_of_gauss_xi,1),
"*",err,error))//
" and it needs to be >= "// &
3407 CALL flagerror(local_error,err,error,*999)
3410 CALL flagerror(
"Quadrature basis is not associated.",err,error,*999)
3413 CALL flagerror(
"Basis has not finished.",err,error,*999)
3416 CALL flagerror(
"Basis is not associated.",err,error,*999)
3419 exits(
"BASIS_QUADRATURE_NUMBER_OF_GAUSS_XI_GET")
3421 999 errorsexits(
"BASIS_QUADRATURE_NUMBER_OF_GAUSS_XI_GET",err,error)
3436 INTEGER(INTG),
INTENT(IN) :: NUMBER_OF_GAUSS_XI(:)
3437 INTEGER(INTG),
INTENT(OUT) :: ERR
3443 enters(
"BASIS_QUADRATURE_NUMBER_OF_GAUSS_XI_SET",err,error,*999)
3445 IF(
ASSOCIATED(basis))
THEN 3446 IF(basis%BASIS_FINISHED)
THEN 3447 CALL flagerror(
"Basis has been finished.",err,error,*999)
3449 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 3450 IF(
SIZE(number_of_gauss_xi,1)==basis%NUMBER_OF_XI)
THEN 3451 IF(any(number_of_gauss_xi<1))
CALL flagerror(
"Invalid number of gauss values.",err,error,*999)
3452 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(1:basis%NUMBER_OF_XI)=number_of_gauss_xi(1:basis%NUMBER_OF_XI)
3454 DO ni=1,basis%NUMBER_OF_XI
3455 SELECT CASE(basis%INTERPOLATION_XI(ni))
3457 IF(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)<2)
THEN 3459 &
" Gauss points are insufficient for linear Lagrange interpolation." 3463 IF(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)<2)
THEN 3465 &
" Gauss points are insufficient for quadratic Lagrange interpolation." 3469 IF(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)<3)
THEN 3471 &
" Gauss points are insufficient for cubic Lagrange interpolation." 3475 IF(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)<2)
THEN 3477 &
" Gauss points are insufficient for quadratic Hermite interpolation." 3481 IF(basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)<3)
THEN 3483 &
" Gauss points are insufficient for cubic Hermite interpolation." 3487 local_warning=
"For simplex elements please set quadrature order rather than number of gauss points." 3490 local_warning=
"For simplex elements please set quadrature order rather than number of gauss points." 3493 local_error=
"Interpolation xi value "//
trim(
number_to_vstring(basis%INTERPOLATION_XI(ni),
"*",err,error))// &
3495 CALL flagerror(local_error,err,error,*999)
3499 local_error=
"The size of the number of Gauss array ("// &
3501 &
") does not match the number of xi directions ("// &
3504 CALL flagerror(local_error,err,error,*999)
3507 CALL flagerror(
"Quadrature basis is not associated.",err,error,*999)
3511 CALL flagerror(
"Basis is not associated.",err,error,*999)
3514 exits(
"BASIS_QUADRATURE_NUMBER_OF_GAUSS_XI_SET")
3516 999 errorsexits(
"BASIS_QUADRATURE_NUMBER_OF_GAUSS_XI_SET",err,error)
3528 INTEGER(INTG),
INTENT(IN) :: SCHEME
3529 INTEGER(INTG),
INTENT(IN) :: GAUSS_POINT
3530 REAL(DP),
INTENT(OUT) :: GAUSS_XI(:)
3531 INTEGER(INTG),
INTENT(OUT) :: ERR
3537 enters(
"BASIS_QUADRATURE_SINGLE_GAUSS_XI_GET",err,error,*999)
3539 IF(
ASSOCIATED(basis))
THEN 3540 IF(basis%BASIS_FINISHED)
THEN 3541 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 3542 quadrature_scheme=>basis%QUADRATURE%QUADRATURE_SCHEME_MAP(scheme)%PTR
3543 IF(
ASSOCIATED(quadrature_scheme))
THEN 3544 IF(
SIZE(gauss_xi)==basis%NUMBER_OF_XI)
THEN 3545 IF(gauss_point>0.AND.gauss_point<=quadrature_scheme%NUMBER_OF_GAUSS)
THEN 3546 gauss_xi(:)=quadrature_scheme%GAUSS_POSITIONS(:,gauss_point)
3548 local_error=
"The specified Gauss point number of "// &
3550 &
"quadrature scheme of the specified element for this field which has "// &
3552 CALL flagerror(local_error,err,error,*999)
3555 local_error=
"The number of xi values to return is invalid and needs to be "// &
3557 CALL flagerror(local_error,err,error,*999)
3560 CALL flagerror(
"The specified quadrature scheme is not associated for this basis.",err,error,*999)
3563 CALL flagerror(
"Quadrature basis is not associated.",err,error,*999)
3566 CALL flagerror(
"Basis has not finished.",err,error,*999)
3569 CALL flagerror(
"Basis is not associated.",err,error,*999)
3572 exits(
"BASIS_QUADRATURE_SINGLE_GAUSS_XI_GET")
3574 999 errorsexits(
"BASIS_QUADRATURE_SINGLE_GAUSS_XI_GET",err,error)
3587 INTEGER(INTG),
INTENT(IN) :: SCHEME
3588 INTEGER(INTG),
INTENT(IN) :: GAUSS_POINTS(:)
3589 REAL(DP),
INTENT(OUT) :: GAUSS_XI(:,:)
3590 INTEGER(INTG),
INTENT(OUT) :: ERR
3593 INTEGER(INTG) :: Gauss_point
3597 enters(
"BASIS_QUADRATURE_MULTIPLE_GAUSS_XI_GET",err,error,*999)
3599 IF(
ASSOCIATED(basis))
THEN 3600 IF(basis%BASIS_FINISHED)
THEN 3601 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 3602 quadrature_scheme=>basis%QUADRATURE%QUADRATURE_SCHEME_MAP(scheme)%PTR
3603 IF(
ASSOCIATED(quadrature_scheme))
THEN 3604 IF(
SIZE(gauss_xi,1)==basis%NUMBER_OF_XI)
THEN 3605 IF(
SIZE(gauss_points)==0)
THEN 3606 IF(
SIZE(gauss_xi,2)==quadrature_scheme%NUMBER_OF_GAUSS)
THEN 3607 gauss_xi=quadrature_scheme%GAUSS_POSITIONS
3609 local_error=
"The number of Gauss Points to return the xi values for is invalid and needs to be "// &
3611 CALL flagerror(local_error,err,error,*999)
3614 DO gauss_point=1,
SIZE(gauss_points)
3615 IF(gauss_points(gauss_point)>0.AND.gauss_points(gauss_point)<=quadrature_scheme%NUMBER_OF_GAUSS)
THEN 3616 gauss_xi(:,gauss_point)=quadrature_scheme%GAUSS_POSITIONS(:,gauss_points(gauss_point))
3618 local_error=
"The specified Gauss point number of "// &
3619 &
trim(
number_to_vstring(gauss_points(gauss_point),
"*",err,error))//
" is invalid for the specified "// &
3620 &
"quadrature scheme of the specified element for this field which has "// &
3622 CALL flagerror(local_error,err,error,*999)
3627 local_error=
"The number of xi values to return is invalid and needs to be "// &
3629 CALL flagerror(local_error,err,error,*999)
3632 CALL flagerror(
"The specified quadrature scheme is not associated for this basis.",err,error,*999)
3635 CALL flagerror(
"Quadrature basis is not associated.",err,error,*999)
3638 CALL flagerror(
"Basis has not finished.",err,error,*999)
3641 CALL flagerror(
"Basis is not associated.",err,error,*999)
3644 exits(
"BASIS_QUADRATURE_MULTIPLE_GAUSS_XI_GET")
3646 999 errorsexits(
"BASIS_QUADRATURE_MULTIPLE_GAUSS_XI_GET",err,error)
3659 INTEGER(INTG),
INTENT(OUT) :: QUADRATURE_ORDER
3660 INTEGER(INTG),
INTENT(OUT) :: ERR
3664 enters(
"BASIS_QUADRATURE_ORDER_GET",err,error,*999)
3666 IF(
ASSOCIATED(basis))
THEN 3667 IF(basis%BASIS_FINISHED)
THEN 3668 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 3669 quadrature_order=basis%QUADRATURE%GAUSS_ORDER
3671 CALL flagerror(
"Quadrature basis is not associated.",err,error,*999)
3674 CALL flagerror(
"Basis has not finished.",err,error,*999)
3677 CALL flagerror(
"Basis is not associated.",err,error,*999)
3680 exits(
"BASIS_QUADRATURE_ORDER_GET")
3682 999 errorsexits(
"BASIS_QUADRATURE_ORDER_GET",err,error)
3694 INTEGER(INTG),
INTENT(IN) :: USER_NUMBER
3695 INTEGER(INTG),
INTENT(IN) :: ORDER
3696 INTEGER(INTG),
INTENT(OUT) :: ERR
3701 enters(
"BASIS_QUADRATURE_ORDER_SET_NUMBER",err,error,*999)
3706 exits(
"BASIS_QUADRATURE_ORDER_SET_NUMBER")
3708 999 errorsexits(
"BASIS_QUADRATURE_ORDER_SET_NUMBER",err,error)
3721 INTEGER(INTG),
INTENT(IN) :: ORDER
3722 INTEGER(INTG),
INTENT(OUT) :: ERR
3727 enters(
"BASIS_QUADRATURE_ORDER_SET_PTR",err,error,*999)
3729 IF(
ASSOCIATED(basis))
THEN 3730 IF(basis%BASIS_FINISHED)
THEN 3731 CALL flagerror(
"Basis has been finished",err,error,*999)
3733 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 3735 IF(order>1.AND.order<=5)
THEN 3736 basis%QUADRATURE%GAUSS_ORDER=order
3739 &
" is invalid. You must specify and order between 1 and 5." 3740 CALL flagerror(local_error,err,error,*999)
3743 CALL flagerror(
"Can only set the quadrature order for simplex basis types.",err,error,*999)
3746 CALL flagerror(
"Quadrature basis is not associated.",err,error,*999)
3750 CALL flagerror(
"Basis is not associated.",err,error,*999)
3753 exits(
"BASIS_QUADRATURE_ORDER_SET_PTR")
3755 999 errorsexits(
"BASIS_QUADRATURE_ORDER_SET_PTR",err,error)
3768 INTEGER(INTG),
INTENT(OUT) :: QUADRATURE_TYPE
3769 INTEGER(INTG),
INTENT(OUT) :: ERR
3773 enters(
"BASIS_QUADRATURE_TYPE_GET",err,error,*999)
3775 IF(
ASSOCIATED(basis))
THEN 3776 IF(basis%BASIS_FINISHED)
THEN 3777 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 3780 CALL flagerror(
"Basis quadrature basis is not associated.",err,error,*999)
3783 CALL flagerror(
"Basis has not finished.",err,error,*999)
3786 CALL flagerror(
"Basis is not associated.",err,error,*999)
3789 exits(
"BASIS_QUADRATURE_TYPE_GET")
3791 999 errorsexits(
"BASIS_QUADRATURE_TYPE_GET",err,error)
3803 INTEGER(INTG),
INTENT(IN) :: USER_NUMBER
3804 INTEGER(INTG),
INTENT(IN) ::
TYPE 3805 INTEGER(INTG),
INTENT(OUT) :: ERR
3810 enters(
"BASIS_QUADRATURE_TYPE_SET_NUMBER",err,error,*999)
3815 exits(
"BASIS_QUADRATURE_TYPE_SET_NUMBER")
3817 999 errorsexits(
"BASIS_QUADRATURE_TYPE_SET_NUMBER",err,error)
3830 INTEGER(INTG),
INTENT(IN) ::
TYPE 3831 INTEGER(INTG),
INTENT(OUT) :: ERR
3836 enters(
"BASIS_QUADRATURE_TYPE_SET_PTR",err,error,*999)
3838 IF(
ASSOCIATED(basis))
THEN 3839 IF(basis%BASIS_FINISHED)
THEN 3840 CALL flagerror(
"Basis has been finished.",err,error,*999)
3842 IF(
ASSOCIATED(basis%QUADRATURE%BASIS))
THEN 3848 CALL flagerror(
"Gauss Laguerre quadrature is not implemented.",err,error,*999)
3851 CALL flagerror(
"Gauss Hermite quadrature is not implemented.",err,error,*999)
3854 CALL flagerror(local_error,err,error,*999)
3857 CALL flagerror(
"Basis quadrature basis is not associated.",err,error,*999)
3861 CALL flagerror(
"Basis is not associated.",err,error,*999)
3864 exits(
"BASIS_QUADRATURE_TYPE_SET_PTR")
3866 999 errorsexits(
"BASIS_QUADRATURE_TYPE_SET_PTR",err,error)
3879 LOGICAL,
INTENT(IN) :: FACE_GAUSS_EVALUATE
3880 INTEGER(INTG),
INTENT(OUT) :: ERR
3883 enters(
"Basis_QuadratureLocalFaceGaussEvaluateSet",err,error,*999)
3885 IF(
ASSOCIATED(basis))
THEN 3886 IF(basis%BASIS_FINISHED)
THEN 3887 CALL flagerror(
"Basis has been finished.",err,error,*999)
3889 basis%QUADRATURE%EVALUATE_FACE_GAUSS=face_gauss_evaluate
3892 CALL flagerror(
"Basis is not associated.",err,error,*999)
3895 exits(
"Basis_QuadratureLocalFaceGaussEvaluateSet")
3897 999 errorsexits(
"Basis_QuadratureLocalFaceGaussEvaluateSet",err,error)
3912 INTEGER(INTG),
INTENT(OUT) :: ERR
3915 INTEGER(INTG) :: MAX_NUM_NODES,ni,nn,ns
3916 INTEGER(INTG),
ALLOCATABLE :: NODES_IN_FACE(:)
3919 enters(
"BASIS_SIMPLEX_BASIS_CREATE",err,error,*999)
3921 IF(
ASSOCIATED(basis))
THEN 3923 basis%NUMBER_OF_XI_COORDINATES=basis%NUMBER_OF_XI+1
3924 ALLOCATE(basis%INTERPOLATION_TYPE(basis%NUMBER_OF_XI_COORDINATES),stat=err)
3925 IF(err/=0)
CALL flagerror(
"Could not allocate INTERPOLATION_TYPE array.",err,error,*999)
3926 ALLOCATE(basis%INTERPOLATION_ORDER(basis%NUMBER_OF_XI_COORDINATES),stat=err)
3927 IF(err/=0)
CALL flagerror(
"Could not allocate INTERPOLATION_ORDER array.",err,error,*999)
3928 ALLOCATE(basis%NUMBER_OF_NODES_XIC(basis%NUMBER_OF_XI_COORDINATES),stat=err)
3929 IF(err/=0)
CALL flagerror(
"Could not allocate NUMBER_OF_NODES_XIC array.",err,error,*999)
3930 basis%DEGENERATE=.false.
3931 basis%NUMBER_OF_COLLAPSED_XI=0
3932 SELECT CASE(basis%NUMBER_OF_XI)
3934 basis%NUMBER_OF_PARTIAL_DERIVATIVES=3
3935 SELECT CASE(basis%INTERPOLATION_XI(1))
3941 basis%NUMBER_OF_NODES_XIC(1)=2
3942 basis%NUMBER_OF_NODES_XIC(2)=2
3944 basis%NUMBER_OF_NODES=2
3950 basis%NUMBER_OF_NODES_XIC(1)=3
3951 basis%NUMBER_OF_NODES_XIC(2)=3
3953 basis%NUMBER_OF_NODES=3
3959 basis%NUMBER_OF_NODES_XIC(1)=4
3960 basis%NUMBER_OF_NODES_XIC(2)=4
3962 basis%NUMBER_OF_NODES=4
3964 CALL flagerror(
"Invalid interpolation type.",err,error,*999)
3967 basis%NUMBER_OF_PARTIAL_DERIVATIVES=6
3968 SELECT CASE(basis%INTERPOLATION_XI(2))
3976 basis%NUMBER_OF_NODES_XIC(1)=2
3977 basis%NUMBER_OF_NODES_XIC(2)=2
3978 basis%NUMBER_OF_NODES_XIC(3)=2
3980 basis%NUMBER_OF_NODES=3
3988 basis%NUMBER_OF_NODES_XIC(1)=3
3989 basis%NUMBER_OF_NODES_XIC(2)=3
3990 basis%NUMBER_OF_NODES_XIC(3)=3
3992 basis%NUMBER_OF_NODES=6
4000 basis%NUMBER_OF_NODES_XIC(1)=4
4001 basis%NUMBER_OF_NODES_XIC(2)=4
4002 basis%NUMBER_OF_NODES_XIC(3)=4
4004 basis%NUMBER_OF_NODES=10
4006 CALL flagerror(
"Invalid interpolation type.",err,error,*999)
4009 basis%NUMBER_OF_PARTIAL_DERIVATIVES=11
4010 SELECT CASE(basis%INTERPOLATION_XI(3))
4020 basis%NUMBER_OF_NODES_XIC(1)=2
4021 basis%NUMBER_OF_NODES_XIC(2)=2
4022 basis%NUMBER_OF_NODES_XIC(3)=2
4023 basis%NUMBER_OF_NODES_XIC(4)=2
4025 basis%NUMBER_OF_NODES=4
4035 basis%NUMBER_OF_NODES_XIC(1)=3
4036 basis%NUMBER_OF_NODES_XIC(2)=3
4037 basis%NUMBER_OF_NODES_XIC(3)=3
4038 basis%NUMBER_OF_NODES_XIC(4)=3
4040 basis%NUMBER_OF_NODES=10
4050 basis%NUMBER_OF_NODES_XIC(1)=4
4051 basis%NUMBER_OF_NODES_XIC(2)=4
4052 basis%NUMBER_OF_NODES_XIC(3)=4
4053 basis%NUMBER_OF_NODES_XIC(4)=4
4055 basis%NUMBER_OF_NODES=20
4057 CALL flagerror(
"Invalid interpolation type.",err,error,*999)
4060 CALL flagerror(
"Invalid number of xi directions.",err,error,*999)
4063 ALLOCATE(basis%NODE_AT_COLLAPSE(basis%NUMBER_OF_NODES),stat=err)
4064 IF(err/=0)
CALL flagerror(
"Could not allocate node at collapse.",err,error,*999)
4065 basis%NODE_AT_COLLAPSE=.false.
4067 ALLOCATE(basis%NODE_POSITION_INDEX(basis%NUMBER_OF_NODES,basis%NUMBER_OF_XI_COORDINATES),stat=err)
4068 IF(err/=0)
CALL flagerror(
"Could not allocate NODE_POSITION_INDEX.",err,error,*999)
4069 SELECT CASE(basis%NUMBER_OF_XI_COORDINATES)
4071 ALLOCATE(basis%NODE_POSITION_INDEX_INV(max_num_nodes,max_num_nodes,1,1),stat=err)
4073 ALLOCATE(basis%NODE_POSITION_INDEX_INV(max_num_nodes,max_num_nodes,max_num_nodes,1),stat=err)
4075 ALLOCATE(basis%NODE_POSITION_INDEX_INV(max_num_nodes,max_num_nodes,max_num_nodes,max_num_nodes),stat=err)
4077 CALL flagerror(
"Invalid number of coordinates.",err,error,*999)
4079 IF(err/=0)
CALL flagerror(
"Could not allocate NODE_POSITION_INDEX_INV.",err,error,*999)
4080 basis%NODE_POSITION_INDEX_INV=0
4083 SELECT CASE(basis%NUMBER_OF_XI)
4085 SELECT CASE(basis%INTERPOLATION_ORDER(1))
4088 basis%NODE_POSITION_INDEX(1,1)=2
4089 basis%NODE_POSITION_INDEX(1,2)=1
4090 basis%NODE_POSITION_INDEX_INV(2,1,1,1)=1
4092 basis%NODE_POSITION_INDEX(2,1)=1
4093 basis%NODE_POSITION_INDEX(2,2)=2
4094 basis%NODE_POSITION_INDEX_INV(1,2,1,1)=2
4097 basis%NODE_POSITION_INDEX(1,1)=3
4098 basis%NODE_POSITION_INDEX(1,2)=1
4099 basis%NODE_POSITION_INDEX_INV(3,1,1,1)=1
4101 basis%NODE_POSITION_INDEX(2,1)=2
4102 basis%NODE_POSITION_INDEX(2,2)=2
4103 basis%NODE_POSITION_INDEX_INV(2,2,1,1)=2
4105 basis%NODE_POSITION_INDEX(3,1)=1
4106 basis%NODE_POSITION_INDEX(3,2)=3
4107 basis%NODE_POSITION_INDEX_INV(1,3,1,1)=3
4110 basis%NODE_POSITION_INDEX(1,1)=4
4111 basis%NODE_POSITION_INDEX(1,2)=1
4112 basis%NODE_POSITION_INDEX_INV(4,1,1,1)=1
4114 basis%NODE_POSITION_INDEX(2,1)=3
4115 basis%NODE_POSITION_INDEX(2,2)=2
4116 basis%NODE_POSITION_INDEX_INV(3,2,1,1)=2
4118 basis%NODE_POSITION_INDEX(3,1)=2
4119 basis%NODE_POSITION_INDEX(3,2)=3
4120 basis%NODE_POSITION_INDEX_INV(2,3,1,1)=3
4122 basis%NODE_POSITION_INDEX(4,1)=1
4123 basis%NODE_POSITION_INDEX(4,2)=4
4124 basis%NODE_POSITION_INDEX_INV(1,4,1,1)=4
4126 CALL flagerror(
"Invalid interpolation order.",err,error,*999)
4129 SELECT CASE(basis%INTERPOLATION_ORDER(1))
4132 basis%NODE_POSITION_INDEX(1,1)=2
4133 basis%NODE_POSITION_INDEX(1,2)=1
4134 basis%NODE_POSITION_INDEX(1,3)=1
4135 basis%NODE_POSITION_INDEX_INV(2,1,1,1)=1
4137 basis%NODE_POSITION_INDEX(2,1)=1
4138 basis%NODE_POSITION_INDEX(2,2)=2
4139 basis%NODE_POSITION_INDEX(2,3)=1
4140 basis%NODE_POSITION_INDEX_INV(1,2,1,1)=2
4142 basis%NODE_POSITION_INDEX(3,1)=1
4143 basis%NODE_POSITION_INDEX(3,2)=1
4144 basis%NODE_POSITION_INDEX(3,3)=2
4145 basis%NODE_POSITION_INDEX_INV(1,1,2,1)=3
4148 basis%NODE_POSITION_INDEX(1,1)=3
4149 basis%NODE_POSITION_INDEX(1,2)=1
4150 basis%NODE_POSITION_INDEX(1,3)=1
4151 basis%NODE_POSITION_INDEX_INV(3,1,1,1)=1
4153 basis%NODE_POSITION_INDEX(2,1)=1
4154 basis%NODE_POSITION_INDEX(2,2)=3
4155 basis%NODE_POSITION_INDEX(2,3)=1
4156 basis%NODE_POSITION_INDEX_INV(1,3,1,1)=2
4158 basis%NODE_POSITION_INDEX(3,1)=1
4159 basis%NODE_POSITION_INDEX(3,2)=1
4160 basis%NODE_POSITION_INDEX(3,3)=3
4161 basis%NODE_POSITION_INDEX_INV(1,1,3,1)=3
4163 basis%NODE_POSITION_INDEX(4,1)=2
4164 basis%NODE_POSITION_INDEX(4,2)=2
4165 basis%NODE_POSITION_INDEX(4,3)=1
4166 basis%NODE_POSITION_INDEX_INV(2,2,1,1)=4
4168 basis%NODE_POSITION_INDEX(5,1)=1
4169 basis%NODE_POSITION_INDEX(5,2)=2
4170 basis%NODE_POSITION_INDEX(5,3)=2
4171 basis%NODE_POSITION_INDEX_INV(1,2,2,1)=5
4173 basis%NODE_POSITION_INDEX(6,1)=2
4174 basis%NODE_POSITION_INDEX(6,2)=1
4175 basis%NODE_POSITION_INDEX(6,3)=2
4176 basis%NODE_POSITION_INDEX_INV(2,1,2,1)=6
4179 basis%NODE_POSITION_INDEX(1,1)=4
4180 basis%NODE_POSITION_INDEX(1,2)=1
4181 basis%NODE_POSITION_INDEX(1,3)=1
4182 basis%NODE_POSITION_INDEX_INV(4,1,1,1)=1
4184 basis%NODE_POSITION_INDEX(2,1)=1
4185 basis%NODE_POSITION_INDEX(2,2)=4
4186 basis%NODE_POSITION_INDEX(2,3)=1
4187 basis%NODE_POSITION_INDEX_INV(1,4,1,1)=2
4189 basis%NODE_POSITION_INDEX(3,1)=1
4190 basis%NODE_POSITION_INDEX(3,2)=1
4191 basis%NODE_POSITION_INDEX(3,3)=4
4192 basis%NODE_POSITION_INDEX_INV(1,1,4,1)=3
4194 basis%NODE_POSITION_INDEX(4,1)=3
4195 basis%NODE_POSITION_INDEX(4,2)=2
4196 basis%NODE_POSITION_INDEX(4,3)=1
4197 basis%NODE_POSITION_INDEX_INV(3,2,1,1)=4
4199 basis%NODE_POSITION_INDEX(5,1)=2
4200 basis%NODE_POSITION_INDEX(5,2)=3
4201 basis%NODE_POSITION_INDEX(5,3)=1
4202 basis%NODE_POSITION_INDEX_INV(2,3,1,1)=5
4204 basis%NODE_POSITION_INDEX(6,1)=1
4205 basis%NODE_POSITION_INDEX(6,2)=3
4206 basis%NODE_POSITION_INDEX(6,3)=2
4207 basis%NODE_POSITION_INDEX_INV(1,3,2,1)=6
4209 basis%NODE_POSITION_INDEX(7,1)=1
4210 basis%NODE_POSITION_INDEX(7,2)=2
4211 basis%NODE_POSITION_INDEX(7,3)=3
4212 basis%NODE_POSITION_INDEX_INV(1,2,3,1)=7
4214 basis%NODE_POSITION_INDEX(8,1)=2
4215 basis%NODE_POSITION_INDEX(8,2)=1
4216 basis%NODE_POSITION_INDEX(8,3)=3
4217 basis%NODE_POSITION_INDEX_INV(2,1,3,1)=8
4219 basis%NODE_POSITION_INDEX(9,1)=3
4220 basis%NODE_POSITION_INDEX(9,2)=1
4221 basis%NODE_POSITION_INDEX(9,3)=2
4222 basis%NODE_POSITION_INDEX_INV(3,1,2,1)=9
4224 basis%NODE_POSITION_INDEX(10,1)=2
4225 basis%NODE_POSITION_INDEX(10,2)=2
4226 basis%NODE_POSITION_INDEX(10,3)=2
4227 basis%NODE_POSITION_INDEX_INV(2,2,2,1)=10
4229 CALL flagerror(
"Invalid interpolation order",err,error,*999)
4232 SELECT CASE(basis%INTERPOLATION_ORDER(1))
4235 basis%NODE_POSITION_INDEX(1,1)=2
4236 basis%NODE_POSITION_INDEX(1,2)=1
4237 basis%NODE_POSITION_INDEX(1,3)=1
4238 basis%NODE_POSITION_INDEX(1,4)=1
4239 basis%NODE_POSITION_INDEX_INV(2,1,1,1)=1
4241 basis%NODE_POSITION_INDEX(2,1)=1
4242 basis%NODE_POSITION_INDEX(2,2)=2
4243 basis%NODE_POSITION_INDEX(2,3)=1
4244 basis%NODE_POSITION_INDEX(2,4)=1
4245 basis%NODE_POSITION_INDEX_INV(1,2,1,1)=2
4247 basis%NODE_POSITION_INDEX(3,1)=1
4248 basis%NODE_POSITION_INDEX(3,2)=1
4249 basis%NODE_POSITION_INDEX(3,3)=2
4250 basis%NODE_POSITION_INDEX(3,4)=1
4251 basis%NODE_POSITION_INDEX_INV(1,1,2,1)=3
4253 basis%NODE_POSITION_INDEX(4,1)=1
4254 basis%NODE_POSITION_INDEX(4,2)=1
4255 basis%NODE_POSITION_INDEX(4,3)=1
4256 basis%NODE_POSITION_INDEX(4,4)=2
4257 basis%NODE_POSITION_INDEX_INV(1,1,1,2)=4
4259 ALLOCATE(nodes_in_face(12),stat=err)
4260 IF(err/=0)
CALL flagerror(
"Could not allocate NODES_IN_FACE.",err,error,*999)
4261 nodes_in_face(:)=[2,3,4,1,3,4,1,2,4,1,2,3]
4265 basis%NODE_POSITION_INDEX(1,1)=3
4266 basis%NODE_POSITION_INDEX(1,2)=1
4267 basis%NODE_POSITION_INDEX(1,3)=1
4268 basis%NODE_POSITION_INDEX(1,4)=1
4269 basis%NODE_POSITION_INDEX_INV(3,1,1,1)=1
4271 basis%NODE_POSITION_INDEX(2,1)=1
4272 basis%NODE_POSITION_INDEX(2,2)=3
4273 basis%NODE_POSITION_INDEX(2,3)=1
4274 basis%NODE_POSITION_INDEX(2,4)=1
4275 basis%NODE_POSITION_INDEX_INV(1,3,1,1)=2
4277 basis%NODE_POSITION_INDEX(3,1)=1
4278 basis%NODE_POSITION_INDEX(3,2)=1
4279 basis%NODE_POSITION_INDEX(3,3)=3
4280 basis%NODE_POSITION_INDEX(3,4)=1
4281 basis%NODE_POSITION_INDEX_INV(1,1,3,1)=3
4283 basis%NODE_POSITION_INDEX(4,1)=1
4284 basis%NODE_POSITION_INDEX(4,2)=1
4285 basis%NODE_POSITION_INDEX(4,3)=1
4286 basis%NODE_POSITION_INDEX(4,4)=3
4287 basis%NODE_POSITION_INDEX_INV(1,1,1,3)=4
4289 basis%NODE_POSITION_INDEX(5,1)=2
4290 basis%NODE_POSITION_INDEX(5,2)=2
4291 basis%NODE_POSITION_INDEX(5,3)=1
4292 basis%NODE_POSITION_INDEX(5,4)=1
4293 basis%NODE_POSITION_INDEX_INV(2,2,1,1)=5
4295 basis%NODE_POSITION_INDEX(6,1)=2
4296 basis%NODE_POSITION_INDEX(6,2)=1
4297 basis%NODE_POSITION_INDEX(6,3)=2
4298 basis%NODE_POSITION_INDEX(6,4)=1
4299 basis%NODE_POSITION_INDEX_INV(2,1,2,1)=6
4301 basis%NODE_POSITION_INDEX(7,1)=2
4302 basis%NODE_POSITION_INDEX(7,2)=1
4303 basis%NODE_POSITION_INDEX(7,3)=1
4304 basis%NODE_POSITION_INDEX(7,4)=2
4305 basis%NODE_POSITION_INDEX_INV(2,1,1,2)=7
4307 basis%NODE_POSITION_INDEX(8,1)=1
4308 basis%NODE_POSITION_INDEX(8,2)=2
4309 basis%NODE_POSITION_INDEX(8,3)=2
4310 basis%NODE_POSITION_INDEX(8,4)=1
4311 basis%NODE_POSITION_INDEX_INV(1,2,2,1)=8
4313 basis%NODE_POSITION_INDEX(9,1)=1
4314 basis%NODE_POSITION_INDEX(9,2)=1
4315 basis%NODE_POSITION_INDEX(9,3)=2
4316 basis%NODE_POSITION_INDEX(9,4)=2
4317 basis%NODE_POSITION_INDEX_INV(1,1,2,2)=9
4319 basis%NODE_POSITION_INDEX(10,1)=1
4320 basis%NODE_POSITION_INDEX(10,2)=2
4321 basis%NODE_POSITION_INDEX(10,3)=1
4322 basis%NODE_POSITION_INDEX(10,4)=2
4323 basis%NODE_POSITION_INDEX_INV(1,2,1,2)=10
4325 ALLOCATE(nodes_in_face(24),stat=err)
4326 IF(err/=0)
CALL flagerror(
"Could not allocate NODES_IN_FACE.",err,error,*999)
4327 nodes_in_face(:)=[2,3,4,8,9,10,1,3,4,6,9,7,1,2,4,5,10,7,1,2,3,5,8,6]
4331 basis%NODE_POSITION_INDEX(1,1)=4
4332 basis%NODE_POSITION_INDEX(1,2)=1
4333 basis%NODE_POSITION_INDEX(1,3)=1
4334 basis%NODE_POSITION_INDEX(1,4)=1
4335 basis%NODE_POSITION_INDEX_INV(4,1,1,1)=1
4337 basis%NODE_POSITION_INDEX(2,1)=1
4338 basis%NODE_POSITION_INDEX(2,2)=4
4339 basis%NODE_POSITION_INDEX(2,3)=1
4340 basis%NODE_POSITION_INDEX(2,4)=1
4341 basis%NODE_POSITION_INDEX_INV(1,4,1,1)=2
4343 basis%NODE_POSITION_INDEX(3,1)=1
4344 basis%NODE_POSITION_INDEX(3,2)=1
4345 basis%NODE_POSITION_INDEX(3,3)=4
4346 basis%NODE_POSITION_INDEX(3,4)=1
4347 basis%NODE_POSITION_INDEX_INV(1,1,4,1)=3
4349 basis%NODE_POSITION_INDEX(4,1)=1
4350 basis%NODE_POSITION_INDEX(4,2)=1
4351 basis%NODE_POSITION_INDEX(4,3)=1
4352 basis%NODE_POSITION_INDEX(4,4)=4
4353 basis%NODE_POSITION_INDEX_INV(1,1,1,4)=4
4355 basis%NODE_POSITION_INDEX(5,1)=3
4356 basis%NODE_POSITION_INDEX(5,2)=2
4357 basis%NODE_POSITION_INDEX(5,3)=1
4358 basis%NODE_POSITION_INDEX(5,4)=1
4359 basis%NODE_POSITION_INDEX_INV(3,2,1,1)=5
4361 basis%NODE_POSITION_INDEX(6,1)=2
4362 basis%NODE_POSITION_INDEX(6,2)=3
4363 basis%NODE_POSITION_INDEX(6,3)=1
4364 basis%NODE_POSITION_INDEX(6,4)=1
4365 basis%NODE_POSITION_INDEX_INV(2,3,1,1)=6
4367 basis%NODE_POSITION_INDEX(7,1)=3
4368 basis%NODE_POSITION_INDEX(7,2)=1
4369 basis%NODE_POSITION_INDEX(7,3)=2
4370 basis%NODE_POSITION_INDEX(7,4)=1
4371 basis%NODE_POSITION_INDEX_INV(3,1,2,1)=7
4373 basis%NODE_POSITION_INDEX(8,1)=2
4374 basis%NODE_POSITION_INDEX(8,2)=1
4375 basis%NODE_POSITION_INDEX(8,3)=3
4376 basis%NODE_POSITION_INDEX(8,4)=1
4377 basis%NODE_POSITION_INDEX_INV(2,1,3,1)=8
4379 basis%NODE_POSITION_INDEX(9,1)=3
4380 basis%NODE_POSITION_INDEX(9,2)=1
4381 basis%NODE_POSITION_INDEX(9,3)=1
4382 basis%NODE_POSITION_INDEX(9,4)=2
4383 basis%NODE_POSITION_INDEX_INV(3,1,1,2)=9
4385 basis%NODE_POSITION_INDEX(10,1)=2
4386 basis%NODE_POSITION_INDEX(10,2)=1
4387 basis%NODE_POSITION_INDEX(10,3)=1
4388 basis%NODE_POSITION_INDEX(10,4)=3
4389 basis%NODE_POSITION_INDEX_INV(2,1,1,3)=10
4391 basis%NODE_POSITION_INDEX(11,1)=1
4392 basis%NODE_POSITION_INDEX(11,2)=3
4393 basis%NODE_POSITION_INDEX(11,3)=2
4394 basis%NODE_POSITION_INDEX(11,4)=1
4395 basis%NODE_POSITION_INDEX_INV(1,3,2,1)=11
4397 basis%NODE_POSITION_INDEX(12,1)=1
4398 basis%NODE_POSITION_INDEX(12,2)=2
4399 basis%NODE_POSITION_INDEX(12,3)=3
4400 basis%NODE_POSITION_INDEX(12,4)=1
4401 basis%NODE_POSITION_INDEX_INV(1,2,3,1)=12
4403 basis%NODE_POSITION_INDEX(13,1)=1
4404 basis%NODE_POSITION_INDEX(13,2)=1
4405 basis%NODE_POSITION_INDEX(13,3)=3
4406 basis%NODE_POSITION_INDEX(13,4)=2
4407 basis%NODE_POSITION_INDEX_INV(1,1,3,2)=13
4409 basis%NODE_POSITION_INDEX(14,1)=1
4410 basis%NODE_POSITION_INDEX(14,2)=1
4411 basis%NODE_POSITION_INDEX(14,3)=2
4412 basis%NODE_POSITION_INDEX(14,4)=3
4413 basis%NODE_POSITION_INDEX_INV(1,1,2,3)=14
4415 basis%NODE_POSITION_INDEX(15,1)=1
4416 basis%NODE_POSITION_INDEX(15,2)=3
4417 basis%NODE_POSITION_INDEX(15,3)=1
4418 basis%NODE_POSITION_INDEX(15,4)=2
4419 basis%NODE_POSITION_INDEX_INV(1,3,1,2)=15
4421 basis%NODE_POSITION_INDEX(16,1)=1
4422 basis%NODE_POSITION_INDEX(16,2)=2
4423 basis%NODE_POSITION_INDEX(16,3)=1
4424 basis%NODE_POSITION_INDEX(16,4)=3
4425 basis%NODE_POSITION_INDEX_INV(1,2,1,3)=16
4427 basis%NODE_POSITION_INDEX(17,1)=2
4428 basis%NODE_POSITION_INDEX(17,2)=2
4429 basis%NODE_POSITION_INDEX(17,3)=2
4430 basis%NODE_POSITION_INDEX(17,4)=1
4431 basis%NODE_POSITION_INDEX_INV(2,2,2,1)=17
4433 basis%NODE_POSITION_INDEX(18,1)=2
4434 basis%NODE_POSITION_INDEX(18,2)=2
4435 basis%NODE_POSITION_INDEX(18,3)=1
4436 basis%NODE_POSITION_INDEX(18,4)=2
4437 basis%NODE_POSITION_INDEX_INV(2,2,1,2)=18
4439 basis%NODE_POSITION_INDEX(19,1)=2
4440 basis%NODE_POSITION_INDEX(19,2)=1
4441 basis%NODE_POSITION_INDEX(19,3)=2
4442 basis%NODE_POSITION_INDEX(19,4)=2
4443 basis%NODE_POSITION_INDEX_INV(2,1,2,2)=19
4445 basis%NODE_POSITION_INDEX(20,1)=1
4446 basis%NODE_POSITION_INDEX(20,2)=2
4447 basis%NODE_POSITION_INDEX(20,3)=2
4448 basis%NODE_POSITION_INDEX(20,4)=2
4449 basis%NODE_POSITION_INDEX_INV(1,2,2,2)=20
4451 ALLOCATE(nodes_in_face(40),stat=err)
4452 IF(err/=0)
CALL flagerror(
"Could not allocate NODES_IN_FACE.",err,error,*999)
4453 nodes_in_face(:)=[2,3,4,11,12,13,14,16,15,20,1,3,4,7,8,13,14,10,9,&
4454 &19,1,2,4,5,6,15,16,10,9,18,1,2,3,5,6,14,12,8,7,17]
4457 CALL flagerror(
"Invalid interpolation order.",err,error,*999)
4460 CALL flagerror(
"Invalid number of xi directions.",err,error,*999)
4463 basis%MAXIMUM_NUMBER_OF_DERIVATIVES=1
4464 basis%NUMBER_OF_ELEMENT_PARAMETERS=basis%NUMBER_OF_NODES
4466 ALLOCATE(basis%NUMBER_OF_DERIVATIVES(basis%NUMBER_OF_NODES),stat=err)
4467 IF(err/=0)
CALL flagerror(
"Could not allocate NUMBER_OF_DERIVATIVES.",err,error,*999)
4468 ALLOCATE(basis%DERIVATIVE_ORDER_INDEX(basis%MAXIMUM_NUMBER_OF_DERIVATIVES,basis%NUMBER_OF_NODES,basis%NUMBER_OF_XI), &
4470 IF(err/=0)
CALL flagerror(
"Could not allocate DERIVATIVE_ORDER_INDEX.",err,error,*999)
4473 IF(err/=0)
CALL flagerror(
"Could not allocate DERIVATIVE_ORDER_INDEX_INV.",err,error,*999)
4474 ALLOCATE(basis%PARTIAL_DERIVATIVE_INDEX(basis%MAXIMUM_NUMBER_OF_DERIVATIVES,basis%NUMBER_OF_NODES),stat=err)
4475 IF(err/=0)
CALL flagerror(
"Could not allocate PARTIAL_DERIVATIVE_INDEX.",err,error,*999)
4476 ALLOCATE(basis%ELEMENT_PARAMETER_INDEX(basis%MAXIMUM_NUMBER_OF_DERIVATIVES,basis%NUMBER_OF_NODES),stat=err)
4477 IF(err/=0)
CALL flagerror(
"Could not allocate ELEMENT_PARAMETER_INDEX.",err,error,*999)
4478 ALLOCATE(basis%ELEMENT_PARAMETER_INDEX_INV(2,basis%NUMBER_OF_ELEMENT_PARAMETERS),stat=err)
4479 IF(err/=0)
CALL flagerror(
"Could not allocate ELEMENT_PARAMETER_INDEX_INV.",err,error,*999)
4482 basis%DERIVATIVE_ORDER_INDEX_INV=0
4483 DO nn=1,basis%NUMBER_OF_NODES
4484 basis%NUMBER_OF_DERIVATIVES(nn)=1
4485 DO ni=1,basis%NUMBER_OF_XI
4486 basis%DERIVATIVE_ORDER_INDEX(1,nn,ni)=1
4489 basis%ELEMENT_PARAMETER_INDEX(1,nn)=ns
4490 basis%ELEMENT_PARAMETER_INDEX_INV(1,ns)=nn
4491 basis%ELEMENT_PARAMETER_INDEX_INV(2,ns)=1
4493 basis%DERIVATIVE_ORDER_INDEX_INV(basis%DERIVATIVE_ORDER_INDEX(1,nn,1),1,1,nn)=1
4497 SELECT CASE(basis%NUMBER_OF_XI)
4499 basis%NUMBER_OF_LOCAL_LINES=1
4500 ALLOCATE(basis%NUMBER_OF_NODES_IN_LOCAL_LINE(basis%NUMBER_OF_LOCAL_LINES),stat=err)
4501 IF(err/=0)
CALL flagerror(
"Could not allocate number of nodes in local line.",err,error,*999)
4502 ALLOCATE(basis%LOCAL_LINE_XI_DIRECTION(basis%NUMBER_OF_LOCAL_LINES),stat=err)
4503 IF(err/=0)
CALL flagerror(
"Could not allocate local line xi direction.",err,error,*999)
4504 basis%LOCAL_LINE_XI_DIRECTION(1)=1
4505 ALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_LINE(basis%NUMBER_OF_NODES_XIC(1),basis%NUMBER_OF_LOCAL_LINES),stat=err)
4506 IF(err/=0)
CALL flagerror(
"Could not allocate node numbers in local line.",err,error,*999)
4507 ALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE(basis%NUMBER_OF_NODES_XIC(1),basis%NUMBER_OF_LOCAL_LINES),stat=err)
4508 IF(err/=0)
CALL flagerror(
"Could not allocate derivative numbers in local line.",err,error,*999)
4510 ALLOCATE(basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE(basis%NUMBER_OF_NODES_XIC(1)**2,basis%NUMBER_OF_LOCAL_LINES) &
4512 IF(err/=0)
CALL flagerror(
"Could not allocate element parameters in local line.",err,error,*999)
4513 basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE=1
4515 SELECT CASE(basis%INTERPOLATION_ORDER(1))
4518 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(1)=2
4519 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,1)=1
4520 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,1)=2
4523 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(1)=3
4524 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,1)=1
4525 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,1)=2
4526 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,1)=3
4529 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(1)=4
4530 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,1)=1
4531 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,1)=2
4532 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,1)=3
4533 basis%NODE_NUMBERS_IN_LOCAL_LINE(4,1)=4
4536 & err,error))//
" is invalid for a simplex basis type." 4537 CALL flagerror(local_error,err,error,*999)
4542 basis%NUMBER_OF_LOCAL_LINES=3
4543 ALLOCATE(basis%NUMBER_OF_NODES_IN_LOCAL_LINE(basis%NUMBER_OF_LOCAL_LINES),stat=err)
4544 IF(err/=0)
CALL flagerror(
"Could not allocate number of nodes in local line.",err,error,*999)
4545 ALLOCATE(basis%LOCAL_LINE_XI_DIRECTION(basis%NUMBER_OF_LOCAL_LINES),stat=err)
4546 IF(err/=0)
CALL flagerror(
"Could not allocate local line xi direction",err,error,*999)
4547 ALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC),basis%NUMBER_OF_LOCAL_LINES),stat=err)
4548 IF(err/=0)
CALL flagerror(
"Could not allocate node numbers in local line.",err,error,*999)
4549 ALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC),basis%NUMBER_OF_LOCAL_LINES),stat=err)
4550 IF(err/=0)
CALL flagerror(
"Could not allocate derivative numbers in local line.",err,error,*999)
4552 ALLOCATE(basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC)**2,basis%NUMBER_OF_LOCAL_LINES) &
4554 IF(err/=0)
CALL flagerror(
"Could not allocate element parameters in local line.",err,error,*999)
4555 basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE=1
4556 ALLOCATE(basis%LOCAL_XI_NORMAL(basis%NUMBER_OF_LOCAL_LINES),stat=err)
4557 IF(err/=0)
CALL flagerror(
"Could not allocate local line normal.",err,error,*999)
4559 SELECT CASE(basis%INTERPOLATION_ORDER(1))
4562 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(1)=2
4563 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,1)=1
4564 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,1)=2
4565 basis%LOCAL_LINE_XI_DIRECTION(1)=1
4566 basis%LOCAL_XI_NORMAL(1)=3
4568 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(2)=2
4569 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,2)=1
4570 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,2)=3
4571 basis%LOCAL_LINE_XI_DIRECTION(2)=2
4572 basis%LOCAL_XI_NORMAL(2)=2
4574 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(3)=2
4575 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,3)=2
4576 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,3)=3
4577 basis%LOCAL_LINE_XI_DIRECTION(3)=3
4578 basis%LOCAL_XI_NORMAL(3)=1
4581 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(1)=3
4582 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,1)=1
4583 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,1)=4
4584 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,1)=2
4585 basis%LOCAL_LINE_XI_DIRECTION(1)=1
4586 basis%LOCAL_XI_NORMAL(1)=3
4588 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(2)=3
4589 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,2)=1
4590 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,2)=6
4591 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,2)=3
4592 basis%LOCAL_LINE_XI_DIRECTION(2)=2
4593 basis%LOCAL_XI_NORMAL(2)=2
4595 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(3)=3
4596 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,3)=2
4597 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,3)=5
4598 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,3)=3
4599 basis%LOCAL_LINE_XI_DIRECTION(3)=3
4600 basis%LOCAL_XI_NORMAL(3)=1
4603 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(1)=4
4604 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,1)=1
4605 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,1)=4
4606 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,1)=5
4607 basis%NODE_NUMBERS_IN_LOCAL_LINE(4,1)=2
4608 basis%LOCAL_LINE_XI_DIRECTION(1)=1
4609 basis%LOCAL_XI_NORMAL(1)=3
4611 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(2)=4
4612 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,2)=1
4613 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,2)=9
4614 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,2)=8
4615 basis%NODE_NUMBERS_IN_LOCAL_LINE(4,2)=3
4616 basis%LOCAL_LINE_XI_DIRECTION(2)=2
4617 basis%LOCAL_XI_NORMAL(2)=2
4619 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(3)=4
4620 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,3)=2
4621 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,3)=6
4622 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,3)=7
4623 basis%NODE_NUMBERS_IN_LOCAL_LINE(4,3)=3
4624 basis%LOCAL_LINE_XI_DIRECTION(3)=3
4625 basis%LOCAL_XI_NORMAL(3)=1
4627 local_error=
"Interpolation order "//
trim(
number_to_vstring(basis%INTERPOLATION_ORDER(1),
"*",err,error))// &
4628 &
" is invalid for a simplex basis type." 4629 CALL flagerror(local_error,err,error,*999)
4632 basis%NUMBER_OF_LOCAL_LINES=6
4633 basis%NUMBER_OF_LOCAL_FACES=4
4635 ALLOCATE(basis%NUMBER_OF_NODES_IN_LOCAL_LINE(basis%NUMBER_OF_LOCAL_LINES),stat=err)
4636 IF(err/=0)
CALL flagerror(
"Could not allocate number of nodes in local line.",err,error,*999)
4637 ALLOCATE(basis%NUMBER_OF_NODES_IN_LOCAL_FACE(basis%NUMBER_OF_LOCAL_FACES),stat=err)
4638 IF(err/=0)
CALL flagerror(
"Could not allocate number of nodes in local face.",err,error,*999)
4640 ALLOCATE(basis%LOCAL_LINE_XI_DIRECTION(basis%NUMBER_OF_LOCAL_LINES),stat=err)
4641 IF(err/=0)
CALL flagerror(
"Could not allocate local line xi direction.",err,error,*999)
4642 ALLOCATE(basis%LOCAL_FACE_XI_DIRECTION(basis%NUMBER_OF_LOCAL_FACES),stat=err)
4643 IF(err/=0)
CALL flagerror(
"Could not allocate local face xi direction.",err,error,*999)
4645 ALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC),basis%NUMBER_OF_LOCAL_LINES),stat=err)
4646 IF(err/=0)
CALL flagerror(
"Could not allocate node numbers in local line.",err,error,*999)
4647 ALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC),basis%NUMBER_OF_LOCAL_LINES),stat=err)
4648 IF(err/=0)
CALL flagerror(
"Could not allocate derivative numbers in local line.",err,error,*999)
4650 ALLOCATE(basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE(maxval(basis%NUMBER_OF_NODES_XIC)**2,basis%NUMBER_OF_LOCAL_LINES), &
4652 IF(err/=0)
CALL flagerror(
"Could not allocate element parameters in local line.",err,error,*999)
4653 basis%ELEMENT_PARAMETERS_IN_LOCAL_LINE=1
4655 SELECT CASE(basis%INTERPOLATION_ORDER(1))
4657 ALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_FACE(3,basis%NUMBER_OF_LOCAL_FACES),stat=err)
4658 IF(err/=0)
CALL flagerror(
"Could not allocate node numbers in local face.",err,error,*999)
4660 ALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE(0:1,3,basis%NUMBER_OF_LOCAL_FACES),stat=err)
4661 IF(err/=0)
CALL flagerror(
"Could not allocate derivative numbers in local face.",err,error,*999)
4663 ALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_FACE(6,basis%NUMBER_OF_LOCAL_FACES),stat=err)
4664 IF(err/=0)
CALL flagerror(
"Could not allocate node numbers in local face.",err,error,*999)
4666 ALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE(0:1,6,basis%NUMBER_OF_LOCAL_FACES),stat=err)
4667 IF(err/=0)
CALL flagerror(
"Could not allocate derivative numbers in local face.",err,error,*999)
4669 ALLOCATE(basis%NODE_NUMBERS_IN_LOCAL_FACE(10,basis%NUMBER_OF_LOCAL_FACES),stat=err)
4670 IF(err/=0)
CALL flagerror(
"Could not allocate node numbers in local face.",err,error,*999)
4672 ALLOCATE(basis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE(0:1,10,basis%NUMBER_OF_LOCAL_FACES),stat=err)
4673 IF(err/=0)
CALL flagerror(
"Could not allocate derivative numbers in local face.",err,error,*999)
4675 local_error=
"Interpolation order "//
trim(
number_to_vstring(basis%INTERPOLATION_ORDER(1),
"*",err,error))// &
4676 &
" is invalid for a simplex basis type." 4677 CALL flagerror(local_error,err,error,*999)
4680 basis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE(0,:,:)=1
4682 ALLOCATE(basis%LOCAL_XI_NORMAL(basis%NUMBER_OF_LOCAL_LINES),stat=err)
4683 IF(err/=0)
CALL flagerror(
"Could not allocate local line normal.",err,error,*999)
4686 SELECT CASE(basis%INTERPOLATION_ORDER(1))
4689 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(1)=2
4690 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,1)=1
4691 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,1)=2
4692 basis%LOCAL_LINE_XI_DIRECTION(1)=1
4694 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(2)=2
4695 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,2)=1
4696 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,2)=3
4697 basis%LOCAL_LINE_XI_DIRECTION(2)=1
4699 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(3)=2
4700 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,3)=1
4701 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,3)=4
4702 basis%LOCAL_LINE_XI_DIRECTION(3)=1
4704 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(4)=2
4705 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,4)=2
4706 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,4)=3
4707 basis%LOCAL_LINE_XI_DIRECTION(4)=2
4709 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(5)=2
4710 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,5)=2
4711 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,5)=4
4712 basis%LOCAL_LINE_XI_DIRECTION(5)=2
4714 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(6)=2
4715 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,6)=3
4716 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,6)=4
4717 basis%LOCAL_LINE_XI_DIRECTION(6)=3
4719 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(1)=3
4720 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,1)=2
4721 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,1)=3
4722 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,1)=4
4723 basis%LOCAL_FACE_XI_DIRECTION(1)=1
4724 basis%LOCAL_XI_NORMAL(1)=1
4726 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(2)=3
4727 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,2)=1
4728 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,2)=4
4729 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,2)=3
4730 basis%LOCAL_FACE_XI_DIRECTION(2)=2
4731 basis%LOCAL_XI_NORMAL(2)=2
4733 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(3)=3
4734 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,3)=1
4735 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,3)=2
4736 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,3)=4
4737 basis%LOCAL_FACE_XI_DIRECTION(3)=3
4738 basis%LOCAL_XI_NORMAL(3)=3
4740 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(4)=3
4741 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,4)=1
4742 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,4)=3
4743 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,4)=2
4744 basis%LOCAL_FACE_XI_DIRECTION(4)=4
4745 basis%LOCAL_XI_NORMAL(4)=4
4748 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(1)=3
4749 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,1)=1
4750 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,1)=5
4751 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,1)=2
4752 basis%LOCAL_LINE_XI_DIRECTION(1)=1
4754 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(2)=3
4755 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,2)=1
4756 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,2)=6
4757 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,2)=3
4758 basis%LOCAL_LINE_XI_DIRECTION(2)=1
4760 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(3)=3
4761 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,3)=1
4762 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,3)=7
4763 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,3)=4
4764 basis%LOCAL_LINE_XI_DIRECTION(3)=1
4766 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(4)=3
4767 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,4)=2
4768 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,4)=8
4769 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,4)=3
4770 basis%LOCAL_LINE_XI_DIRECTION(4)=2
4772 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(5)=3
4773 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,5)=2
4774 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,5)=10
4775 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,5)=4
4776 basis%LOCAL_LINE_XI_DIRECTION(5)=2
4778 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(6)=3
4779 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,6)=3
4780 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,6)=9
4781 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,6)=4
4782 basis%LOCAL_LINE_XI_DIRECTION(6)=3
4784 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(1)=6
4785 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,1)=2
4786 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,1)=3
4787 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,1)=4
4788 basis%NODE_NUMBERS_IN_LOCAL_FACE(4,1)=8
4789 basis%NODE_NUMBERS_IN_LOCAL_FACE(5,1)=9
4790 basis%NODE_NUMBERS_IN_LOCAL_FACE(6,1)=10
4791 basis%LOCAL_FACE_XI_DIRECTION(1)=1
4792 basis%LOCAL_XI_NORMAL(1)=1
4794 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(2)=6
4795 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,2)=1
4796 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,2)=4
4797 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,2)=3
4798 basis%NODE_NUMBERS_IN_LOCAL_FACE(4,2)=7
4799 basis%NODE_NUMBERS_IN_LOCAL_FACE(5,2)=9
4800 basis%NODE_NUMBERS_IN_LOCAL_FACE(6,2)=6
4801 basis%LOCAL_FACE_XI_DIRECTION(2)=2
4802 basis%LOCAL_XI_NORMAL(2)=2
4804 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(3)=6
4805 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,3)=1
4806 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,3)=2
4807 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,3)=4
4808 basis%NODE_NUMBERS_IN_LOCAL_FACE(4,3)=5
4809 basis%NODE_NUMBERS_IN_LOCAL_FACE(5,3)=10
4810 basis%NODE_NUMBERS_IN_LOCAL_FACE(6,3)=7
4811 basis%LOCAL_FACE_XI_DIRECTION(3)=3
4812 basis%LOCAL_XI_NORMAL(3)=3
4814 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(4)=6
4815 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,4)=1
4816 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,4)=3
4817 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,4)=2
4818 basis%NODE_NUMBERS_IN_LOCAL_FACE(4,4)=6
4819 basis%NODE_NUMBERS_IN_LOCAL_FACE(5,4)=8
4820 basis%NODE_NUMBERS_IN_LOCAL_FACE(6,4)=5
4821 basis%LOCAL_FACE_XI_DIRECTION(4)=4
4822 basis%LOCAL_XI_NORMAL(4)=4
4825 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(1)=4
4826 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,1)=1
4827 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,1)=5
4828 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,1)=6
4829 basis%NODE_NUMBERS_IN_LOCAL_LINE(4,1)=2
4830 basis%LOCAL_LINE_XI_DIRECTION(1)=1
4832 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(2)=4
4833 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,2)=1
4834 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,2)=7
4835 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,2)=8
4836 basis%NODE_NUMBERS_IN_LOCAL_LINE(4,2)=3
4837 basis%LOCAL_LINE_XI_DIRECTION(2)=1
4839 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(3)=4
4840 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,3)=1
4841 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,3)=9
4842 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,3)=10
4843 basis%NODE_NUMBERS_IN_LOCAL_LINE(4,3)=4
4844 basis%LOCAL_LINE_XI_DIRECTION(3)=1
4846 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(4)=4
4847 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,4)=2
4848 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,4)=11
4849 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,4)=12
4850 basis%NODE_NUMBERS_IN_LOCAL_LINE(4,4)=3
4851 basis%LOCAL_LINE_XI_DIRECTION(4)=2
4853 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(5)=4
4854 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,5)=2
4855 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,5)=15
4856 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,5)=16
4857 basis%NODE_NUMBERS_IN_LOCAL_LINE(4,5)=4
4858 basis%LOCAL_LINE_XI_DIRECTION(5)=2
4860 basis%NUMBER_OF_NODES_IN_LOCAL_LINE(6)=4
4861 basis%NODE_NUMBERS_IN_LOCAL_LINE(1,6)=3
4862 basis%NODE_NUMBERS_IN_LOCAL_LINE(2,6)=13
4863 basis%NODE_NUMBERS_IN_LOCAL_LINE(3,6)=14
4864 basis%NODE_NUMBERS_IN_LOCAL_LINE(4,6)=4
4865 basis%LOCAL_LINE_XI_DIRECTION(6)=3
4867 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(1)=10
4868 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,1)=2
4869 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,1)=3
4870 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,1)=4
4871 basis%NODE_NUMBERS_IN_LOCAL_FACE(4,1)=11
4872 basis%NODE_NUMBERS_IN_LOCAL_FACE(5,1)=12
4873 basis%NODE_NUMBERS_IN_LOCAL_FACE(6,1)=13
4874 basis%NODE_NUMBERS_IN_LOCAL_FACE(7,1)=14
4875 basis%NODE_NUMBERS_IN_LOCAL_FACE(8,1)=16
4876 basis%NODE_NUMBERS_IN_LOCAL_FACE(9,1)=15
4877 basis%NODE_NUMBERS_IN_LOCAL_FACE(10,1)=20
4878 basis%LOCAL_FACE_XI_DIRECTION(1)=1
4879 basis%LOCAL_XI_NORMAL(1)=1
4881 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(2)=10
4882 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,2)=1
4883 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,2)=4
4884 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,2)=3
4885 basis%NODE_NUMBERS_IN_LOCAL_FACE(4,2)=9
4886 basis%NODE_NUMBERS_IN_LOCAL_FACE(5,2)=10
4887 basis%NODE_NUMBERS_IN_LOCAL_FACE(6,2)=14
4888 basis%NODE_NUMBERS_IN_LOCAL_FACE(7,2)=13
4889 basis%NODE_NUMBERS_IN_LOCAL_FACE(8,2)=8
4890 basis%NODE_NUMBERS_IN_LOCAL_FACE(9,2)=7
4891 basis%NODE_NUMBERS_IN_LOCAL_FACE(10,2)=19
4892 basis%LOCAL_FACE_XI_DIRECTION(2)=2
4893 basis%LOCAL_XI_NORMAL(2)=2
4895 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(3)=10
4896 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,3)=1
4897 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,3)=2
4898 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,3)=4
4899 basis%NODE_NUMBERS_IN_LOCAL_FACE(4,3)=5
4900 basis%NODE_NUMBERS_IN_LOCAL_FACE(5,3)=6
4901 basis%NODE_NUMBERS_IN_LOCAL_FACE(6,3)=15
4902 basis%NODE_NUMBERS_IN_LOCAL_FACE(7,3)=16
4903 basis%NODE_NUMBERS_IN_LOCAL_FACE(8,3)=10
4904 basis%NODE_NUMBERS_IN_LOCAL_FACE(9,3)=9
4905 basis%NODE_NUMBERS_IN_LOCAL_FACE(10,3)=18
4906 basis%LOCAL_FACE_XI_DIRECTION(3)=3
4907 basis%LOCAL_XI_NORMAL(3)=3
4909 basis%NUMBER_OF_NODES_IN_LOCAL_FACE(4)=10
4910 basis%NODE_NUMBERS_IN_LOCAL_FACE(1,4)=1
4911 basis%NODE_NUMBERS_IN_LOCAL_FACE(2,4)=3
4912 basis%NODE_NUMBERS_IN_LOCAL_FACE(3,4)=2
4913 basis%NODE_NUMBERS_IN_LOCAL_FACE(4,4)=7
4914 basis%NODE_NUMBERS_IN_LOCAL_FACE(5,4)=8
4915 basis%NODE_NUMBERS_IN_LOCAL_FACE(6,4)=12
4916 basis%NODE_NUMBERS_IN_LOCAL_FACE(7,4)=11
4917 basis%NODE_NUMBERS_IN_LOCAL_FACE(8,4)=6
4918 basis%NODE_NUMBERS_IN_LOCAL_FACE(9,4)=5
4919 basis%NODE_NUMBERS_IN_LOCAL_FACE(10,4)=17
4920 basis%LOCAL_FACE_XI_DIRECTION(4)=4
4921 basis%LOCAL_XI_NORMAL(4)=4
4923 local_error=
"Interpolation order "//
trim(
number_to_vstring(basis%INTERPOLATION_ORDER(1),
"*",err,error))// &
4924 &
" is invalid for a simplex basis type." 4925 CALL flagerror(local_error,err,error,*999)
4928 CALL flagerror(
"Invalid number of xi directions.",err,error,*999)
4934 CALL flagerror(
"Basis is not a simplex basis.",err,error,*999)
4937 CALL flagerror(
"Basis is not associated.",err,error,*999)
4940 exits(
"BASIS_SIMPLEX_BASIS_CREATE")
4942 999 errorsexits(
"BASIS_SIMPLEX_BASIS_CREATE",err,error)
4956 INTEGER(INTG),
INTENT(OUT) :: ERR
4959 INTEGER(INTG) :: DUMMY_ERR,ni,ni2,FACE_XI(2),FACE_XI2(2)
4960 LOGICAL :: LINE_BASIS_DONE,FACE_BASIS_DONE
4964 NULLIFY(new_sub_basis)
4966 enters(
"BASIS_SIMPLEX_FAMILY_CREATE",err,error,*999)
4968 IF(
ASSOCIATED(basis))
THEN 4971 IF(basis%NUMBER_OF_XI>1)
THEN 4974 IF (basis%NUMBER_OF_XI .EQ. 2)
THEN 4975 ALLOCATE(basis%LINE_BASES(basis%NUMBER_OF_XI+1),stat=err)
4977 ALLOCATE(basis%LINE_BASES(basis%NUMBER_OF_XI),stat=err)
4980 IF(err/=0)
CALL flagerror(
"Could not allocate basis line bases.",err,error,*999)
4981 DO ni=1,basis%NUMBER_OF_XI
4982 line_basis_done=.false.
4983 NULLIFY(new_sub_basis)
4985 IF(basis%INTERPOLATION_XI(ni2)==basis%INTERPOLATION_XI(ni).AND. &
4986 basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni2)==basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni))
THEN 4987 line_basis_done=.true.
4991 IF(line_basis_done)
THEN 4992 basis%LINE_BASES(ni)%PTR=>basis%LINE_BASES(ni2)%PTR
4998 basis%LINE_BASES(ni)%PTR=>new_sub_basis
5002 IF (basis%NUMBER_OF_XI .EQ. 2)
THEN 5003 basis%LINE_BASES(basis%NUMBER_OF_XI+1)%PTR=>basis%LINE_BASES(basis%NUMBER_OF_XI)%PTR
5006 IF(basis%NUMBER_OF_XI>2)
THEN 5008 ALLOCATE(basis%FACE_BASES(basis%NUMBER_OF_XI),stat=err)
5009 IF(err/=0)
CALL flagerror(
"Could not allocate basis face bases.",err,error,*999)
5010 DO ni=1,basis%NUMBER_OF_XI
5014 face_basis_done=.false.
5015 NULLIFY(new_sub_basis)
5019 IF(basis%INTERPOLATION_XI(face_xi2(1))==basis%INTERPOLATION_XI(face_xi(1)).AND. &
5020 & basis%INTERPOLATION_XI(face_xi2(2))==basis%INTERPOLATION_XI(face_xi(2)).AND. &
5021 & basis%QUADRATURE%NUMBER_OF_GAUSS_XI(face_xi2(1))==basis%QUADRATURE%NUMBER_OF_GAUSS_XI(face_xi(1)).AND. &
5022 & basis%QUADRATURE%NUMBER_OF_GAUSS_XI(face_xi2(2))==basis%QUADRATURE%NUMBER_OF_GAUSS_XI(face_xi(1)))
THEN 5023 face_basis_done=.true.
5027 IF(face_basis_done)
THEN 5028 basis%FACE_BASES(ni)%PTR=>basis%FACE_BASES(ni2)%PTR
5034 new_sub_basis%LINE_BASES(1)%PTR=>basis%LINE_BASES(face_xi(1))%PTR
5035 new_sub_basis%LINE_BASES(2)%PTR=>basis%LINE_BASES(face_xi(2))%PTR
5036 basis%FACE_BASES(ni)%PTR=>new_sub_basis
5040 ALLOCATE(basis%FACE_BASES(1),stat=err)
5041 IF(err/=0)
CALL flagerror(
"Could not allocate basis face bases.",err,error,*999)
5042 basis%FACE_BASES(1)%PTR=>basis
5045 ALLOCATE(basis%LINE_BASES(1),stat=err)
5046 IF(err/=0)
CALL flagerror(
"Could not allocate basis line bases.",err,error,*999)
5047 basis%LINE_BASES(1)%PTR=>basis
5048 NULLIFY(basis%FACE_BASES)
5051 CALL flagerror(
"Basis is not associated.",err,error,*999)
5054 exits(
"BASIS_SIMPLEX_FAMILY_CREATE")
5056 999
IF(
ASSOCIATED(new_sub_basis))
CALL basis_family_destroy(new_sub_basis%USER_NUMBER,new_sub_basis%FAMILY_NUMBER, &
5057 & dummy_err,dummy_error,*998)
5058 998 errorsexits(
"BASIS_SIMPLEX_FAMILY_CREATE",err,error)
5116 INTEGER(INTG),
INTENT(IN) :: NODE_NUMBER
5117 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIV_INDEX
5118 REAL(DP),
INTENT(IN) :: XL(:)
5119 INTEGER(INTG),
INTENT(OUT) :: ERR
5122 REAL(DP) :: BASIS_SIMPLEX_BASIS_EVALUATE
5126 enters(
"BASIS_SIMPLEX_BASIS_EVALUATE",err,error,*999)
5128 basis_simplex_basis_evaluate=0.0_dp
5129 IF(
ASSOCIATED(basis))
THEN 5131 SELECT CASE(basis%NUMBER_OF_XI)
5133 SELECT CASE(partial_deriv_index)
5135 basis_simplex_basis_evaluate= &
5138 basis_simplex_basis_evaluate= &
5141 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5144 basis_simplex_basis_evaluate= &
5147 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5150 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5153 local_error=
"The specified partial derivative index of "//
trim(
number_to_vstring(partial_deriv_index,
"*",err,error))// &
5154 &
" is invalid for a Simplex line basis." 5155 CALL flagerror(local_error,err,error,*999)
5158 SELECT CASE(partial_deriv_index)
5160 basis_simplex_basis_evaluate= &
5163 basis_simplex_basis_evaluate= &
5166 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5169 basis_simplex_basis_evaluate= &
5172 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5175 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5178 basis_simplex_basis_evaluate= &
5181 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5184 basis_simplex_basis_evaluate= &
5187 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5190 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5193 basis_simplex_basis_evaluate= &
5196 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5199 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5202 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5205 local_error=
"The specified partial derivative index of "//
trim(
number_to_vstring(partial_deriv_index,
"*",err,error))// &
5206 &
" is invalid for a Simplex triangle basis." 5207 CALL flagerror(local_error,err,error,*999)
5210 SELECT CASE(partial_deriv_index)
5212 basis_simplex_basis_evaluate= &
5215 basis_simplex_basis_evaluate= &
5218 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5221 basis_simplex_basis_evaluate= &
5224 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5227 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5230 basis_simplex_basis_evaluate= &
5233 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5236 basis_simplex_basis_evaluate= &
5239 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5242 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5245 basis_simplex_basis_evaluate= &
5248 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5251 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5254 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5257 basis_simplex_basis_evaluate= &
5260 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5263 basis_simplex_basis_evaluate= &
5266 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5269 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5272 basis_simplex_basis_evaluate= &
5275 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5278 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5281 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5284 basis_simplex_basis_evaluate= &
5287 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5290 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5293 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5296 basis_simplex_basis_evaluate= &
5299 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5302 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5305 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5308 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5311 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5314 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate+ &
5317 basis_simplex_basis_evaluate=basis_simplex_basis_evaluate- &
5320 local_error=
"The specified partial derivative index of "//
trim(
number_to_vstring(partial_deriv_index,
"*",err,error))// &
5321 &
" is invalid for a Simplex tetrahedra basis." 5322 CALL flagerror(local_error,err,error,*999)
5325 local_error=
"Invalid number of Xi coordinates. The number of xi coordinates for this basis is "// &
5327 CALL flagerror(local_error,err,error,*999)
5331 CALL flagerror(
"Basis is not a simplex basis.",err,error,*999)
5334 CALL flagerror(
"Basis is not associated.",err,error,*999)
5337 exits(
"BASIS_SIMPLEX_BASIS_EVALUATE.")
5339 999 errorsexits(
"BASIS_SIMPLEX_BASIS_EVALUATE.",err,error)
5352 INTEGER(INTG),
INTENT(IN) :: NODE_NUMBER
5353 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIV_INDEX
5354 REAL(DP),
INTENT(IN) :: XL(:)
5355 INTEGER(INTG),
INTENT(OUT) :: ERR
5358 REAL(DP) :: BASIS_SIMPLEX_BASIS_DERIVATIVE_EVALUATE
5360 INTEGER(INTG) :: nic
5363 enters(
"BASIS_SIMPLEX_BASIS_DERIVATIVE_EVALUATE",err,error,*999)
5365 basis_simplex_basis_derivative_evaluate=1.0_dp
5366 IF(
ASSOCIATED(basis))
THEN 5367 DO nic=1,basis%NUMBER_OF_XI_COORDINATES
5368 SELECT CASE(basis%INTERPOLATION_ORDER(nic))
5370 basis_simplex_basis_derivative_evaluate=basis_simplex_basis_derivative_evaluate* &
5374 basis_simplex_basis_derivative_evaluate=basis_simplex_basis_derivative_evaluate* &
5378 basis_simplex_basis_derivative_evaluate=basis_simplex_basis_derivative_evaluate* &
5382 local_error=
"Interpolation order value "//
trim(
number_to_vstring(basis%INTERPOLATION_ORDER(nic),
"*",err,error))// &
5384 CALL flagerror(local_error,err,error,*999)
5389 CALL flagerror(
"Basis is not associated.",err,error,*999)
5392 exits(
"BASIS_SIMPLEX_BASIS_DERIVATIVE_EVALUATE")
5394 999 errorsexits(
"BASIS_SIMPLEX_BASIS_DERIVATIVE_EVALUATE",err,error)
5407 INTEGER(INTG),
INTENT(IN) :: NUMBER_OF_XI
5408 INTEGER(INTG),
INTENT(IN) :: XI_DIRECTIONS(:)
5410 INTEGER(INTG),
INTENT(OUT) :: ERR
5413 INTEGER(INTG) :: basis_idx,ni,NUMBER_COLLAPSED,NUMBER_END_COLLAPSED
5418 NULLIFY(new_sub_basis)
5419 NULLIFY(new_sub_bases)
5421 enters(
"BASIS_SUB_BASIS_CREATE",err,error,*999)
5423 IF(
ASSOCIATED(parent_basis))
THEN 5424 IF(
ASSOCIATED(sub_basis))
THEN 5425 CALL flagerror(
"The sub-basis is already associated",err,error,*999)
5427 IF(number_of_xi>0.AND.number_of_xi<4)
THEN 5428 IF(any(xi_directions<1).OR.any(xi_directions>3))
CALL flagerror(
"Invalid xi directions specified",err,error,*999)
5429 IF(
SIZE(xi_directions,1)/=number_of_xi) &
5430 &
CALL flagerror(
"The size of the xi directions array must be the same as the number of xi directions",err,error,*999)
5431 ALLOCATE(new_sub_basis,stat=err)
5432 IF(err/=0)
CALL flagerror(
"Could not allocate sub-basis",err,error,*999)
5433 new_sub_basis%USER_NUMBER=parent_basis%USER_NUMBER
5434 new_sub_basis%GLOBAL_NUMBER=parent_basis%GLOBAL_NUMBER
5435 new_sub_basis%FAMILY_NUMBER=parent_basis%NUMBER_OF_SUB_BASES+1
5436 new_sub_basis%NUMBER_OF_SUB_BASES=0
5437 NULLIFY(new_sub_basis%SUB_BASES)
5438 new_sub_basis%PARENT_BASIS=>parent_basis
5439 new_sub_basis%NUMBER_OF_XI=number_of_xi
5440 new_sub_basis%TYPE=parent_basis%TYPE
5441 ALLOCATE(new_sub_basis%INTERPOLATION_XI(number_of_xi),stat=err)
5442 IF(err/=0)
CALL flagerror(
"Could not allocate sub-basis interpolation xi",err,error,*999)
5443 ALLOCATE(new_sub_basis%COLLAPSED_XI(number_of_xi),stat=err)
5444 IF(err/=0)
CALL flagerror(
"Could not allocate sub-basis collapsed xi",err,error,*999)
5446 number_end_collapsed=0
5447 DO ni=1,number_of_xi
5448 new_sub_basis%INTERPOLATION_XI(ni)=parent_basis%INTERPOLATION_XI(xi_directions(ni))
5449 new_sub_basis%COLLAPSED_XI(ni)=parent_basis%COLLAPSED_XI(xi_directions(ni))
5451 number_collapsed=number_collapsed+1
5454 number_end_collapsed=number_end_collapsed+1
5457 IF(number_collapsed==0.OR.number_end_collapsed==0) new_sub_basis%COLLAPSED_XI(1:number_of_xi)=
basis_not_collapsed 5458 NULLIFY(new_sub_basis%QUADRATURE%BASIS)
5460 new_sub_basis%QUADRATURE%TYPE=parent_basis%QUADRATURE%TYPE
5461 DO ni=1,number_of_xi
5462 new_sub_basis%QUADRATURE%NUMBER_OF_GAUSS_XI(ni)=parent_basis%QUADRATURE%NUMBER_OF_GAUSS_XI(xi_directions(ni))
5464 new_sub_basis%BASIS_FINISHED=.true.
5465 IF(number_of_xi>1)
THEN 5466 ALLOCATE(new_sub_basis%LINE_BASES(number_of_xi),stat=err)
5467 IF(err/=0)
CALL flagerror(
"Could not allocate sub-basis line bases",err,error,*999)
5468 IF(number_of_xi>2)
THEN 5469 ALLOCATE(new_sub_basis%FACE_BASES(number_of_xi),stat=err)
5470 IF(err/=0)
CALL flagerror(
"Could not allocate sub-basis face bases",err,error,*999)
5472 ALLOCATE(new_sub_basis%FACE_BASES(1),stat=err)
5473 IF(err/=0)
CALL flagerror(
"Could not allocate sub-basis face bases",err,error,*999)
5474 new_sub_basis%FACE_BASES(1)%PTR=>new_sub_basis
5477 ALLOCATE(new_sub_basis%LINE_BASES(1),stat=err)
5478 IF(err/=0)
CALL flagerror(
"Could not allocate basis line bases",err,error,*999)
5479 new_sub_basis%LINE_BASES(1)%PTR=>new_sub_basis
5480 NULLIFY(new_sub_basis%FACE_BASES)
5483 ALLOCATE(new_sub_bases(parent_basis%NUMBER_OF_SUB_BASES+1),stat=err)
5484 IF(err/=0)
CALL flagerror(
"Could not allocate new sub-bases",err,error,*999)
5485 DO basis_idx=1,parent_basis%NUMBER_OF_SUB_BASES
5486 new_sub_bases(basis_idx)%PTR=>parent_basis%SUB_BASES(basis_idx)%PTR
5488 new_sub_bases(parent_basis%NUMBER_OF_SUB_BASES+1)%PTR=>new_sub_basis
5489 parent_basis%NUMBER_OF_SUB_BASES=parent_basis%NUMBER_OF_SUB_BASES+1
5490 IF(
ASSOCIATED(parent_basis%SUB_BASES))
DEALLOCATE(parent_basis%SUB_BASES)
5491 parent_basis%SUB_BASES=>new_sub_bases
5492 sub_basis=>new_sub_basis
5494 local_error=
"Invalid number of xi directions specified ("// &
5495 &
trim(
number_to_vstring(number_of_xi,
"*",err,error))//
"). You must specify between 1 and 3 xi directions" 5496 CALL flagerror(local_error,err,error,*999)
5500 CALL flagerror(
"Parent basis is not associated",err,error,*999)
5503 exits(
"BASIS_SUB_BASIS_CREATE")
5505 999 errorsexits(
"BASIS_SUB_BASIS_CREATE",err,error)
5518 INTEGER(INTG),
INTENT(OUT) ::
TYPE 5519 INTEGER(INTG),
INTENT(OUT) :: ERR
5523 enters(
"BASIS_TYPE_GET",err,error,*999)
5525 IF(
ASSOCIATED(basis))
THEN 5526 IF(basis%BASIS_FINISHED)
THEN 5529 CALL flagerror(
"Basis has not been finished yet",err,error,*999)
5532 CALL flagerror(
"Basis is not associated",err,error,*999)
5535 exits(
"BASIS_TYPE_GET")
5537 999 errorsexits(
"BASIS_TYPE_GET",err,error)
5549 INTEGER(INTG),
INTENT(IN) :: USER_NUMBER
5550 INTEGER(INTG),
INTENT(IN) ::
TYPE 5551 INTEGER(INTG),
INTENT(OUT) :: ERR
5556 enters(
"BASIS_TYPE_SET_NUMBER",err,error,*999)
5561 exits(
"BASIS_TYPE_SET_NUMBER")
5563 999 errorsexits(
"BASIS_TYPE_SET_NUMBER",err,error)
5576 INTEGER(INTG),
INTENT(IN) ::
TYPE 5577 INTEGER(INTG),
INTENT(OUT) :: ERR
5582 enters(
"BASIS_TYPE_SET_PTR",err,error,*999)
5584 IF(
ASSOCIATED(basis))
THEN 5585 IF(basis%BASIS_FINISHED)
THEN 5586 CALL flagerror(
"Basis has been finished",err,error,*999)
5597 NULLIFY(basis%QUADRATURE%BASIS)
5600 local_error=
"Basis type "//
trim(
number_to_vstring(
TYPE,
"*",ERR,ERROR))//
" is invalid or not implemented" 5601 CALL flagerror(local_error,err,error,*999)
5605 CALL flagerror(
"Basis is not associated",err,error,*999)
5608 exits(
"BASIS_TYPE_SET_PTR")
5610 999 errorsexits(
"BASIS_TYPE_SET_PTR",err,error)
5623 INTEGER(INTG),
INTENT(OUT) :: COLLAPSED_XI(:)
5624 INTEGER(INTG),
INTENT(OUT) :: ERR
5629 enters(
"BASIS_COLLAPSED_XI_GET",err,error,*999)
5631 IF(
ASSOCIATED(basis))
THEN 5632 IF(basis%BASIS_FINISHED)
THEN 5633 IF(
SIZE(collapsed_xi,1)>=
SIZE(basis%COLLAPSED_XI))
THEN 5634 collapsed_xi=basis%COLLAPSED_XI
5636 local_error=
"The size of COLLAPSED_XI is too small. The supplied size is "// &
5639 CALL flagerror(local_error,err,error,*999)
5642 CALL flagerror(
"Basis has not been finished.",err,error,*999)
5645 CALL flagerror(
"Basis is not associated.",err,error,*999)
5648 exits(
"BASIS_COLLAPSED_XI_GET")
5650 999 errorsexits(
"BASIS_COLLAPSED_XI_GET",err,error)
5662 INTEGER(INTG),
INTENT(IN) :: USER_NUMBER
5663 INTEGER(INTG),
INTENT(IN) :: COLLAPSED_XI(:)
5664 INTEGER(INTG),
INTENT(OUT) :: ERR
5669 enters(
"BASIS_COLLAPSED_XI_SET_NUMBER",err,error,*999)
5674 exits(
"BASIS_COLLAPSED_XI_SET_NUMBER")
5676 999 errorsexits(
"BASIS_COLLAPSED_XI_SET_NUMBER",err,error)
5689 INTEGER(INTG),
INTENT(IN) :: COLLAPSED_XI(:)
5690 INTEGER(INTG),
INTENT(OUT) :: ERR
5693 INTEGER(INTG) :: ni1,ni2,ni3,NUMBER_COLLAPSED,COLLAPSED_XI_DIR(3)
5696 enters(
"BASIS_COLLAPSED_XI_SET_PTR",err,error,*999)
5698 IF(
ASSOCIATED(basis))
THEN 5699 IF(basis%BASIS_FINISHED)
THEN 5700 CALL flagerror(
"Basis has been finished",err,error,*999)
5703 IF(basis%NUMBER_OF_XI>1)
THEN 5704 IF(
SIZE(collapsed_xi,1)==basis%NUMBER_OF_XI)
THEN 5706 DO ni1=1,basis%NUMBER_OF_XI
5707 SELECT CASE(collapsed_xi(ni1))
5709 number_collapsed=number_collapsed+1
5710 collapsed_xi_dir(number_collapsed)=ni1
5716 CALL flagerror(local_error,err,error,*999)
5719 IF(number_collapsed>0)
THEN 5720 IF(number_collapsed<basis%NUMBER_OF_XI)
THEN 5721 IF(basis%NUMBER_OF_XI==2)
THEN 5723 ni1=collapsed_xi_dir(1)
5732 local_error=
"Invalid collapsing of a two dimensional basis. Xi direction "// &
5735 CALL flagerror(local_error,err,error,*999)
5739 IF(number_collapsed==1)
THEN 5741 ni1=collapsed_xi_dir(1)
5752 local_error=
"Invalid collapsing of a three dimensional basis. Xi direction "// &
5756 CALL flagerror(local_error,err,error,*999)
5766 local_error=
"Invalid collapsing of a three dimensional basis. Xi direction "// &
5770 CALL flagerror(local_error,err,error,*999)
5773 local_error=
"Invalid collapsing of a three dimensional basis. Xi direction "// &
5777 CALL flagerror(local_error,err,error,*999)
5781 ni1=collapsed_xi_dir(1)
5782 ni2=collapsed_xi_dir(2)
5791 local_error=
"Invalid collapsing of a three dimensional basis. Xi directions "// &
5795 CALL flagerror(local_error,err,error,*999)
5800 local_error=
"Invalid collapsing of basis. The number of collapsed directions ("// &
5802 &
") must be less than the number of xi directions ("// &
5804 CALL flagerror(local_error,err,error,*999)
5808 DO ni1=1,basis%NUMBER_OF_XI
5815 basis%COLLAPSED_XI(1:basis%NUMBER_OF_XI)=collapsed_xi(1:basis%NUMBER_OF_XI)
5817 local_error=
"The size of the xi collapsed array ("// &
5818 &
trim(
number_to_vstring(
SIZE(collapsed_xi,1),
"*",err,error))//
") does not match the number of xi directions ("// &
5821 CALL flagerror(local_error,err,error,*999)
5824 CALL flagerror(
"Can not collapse a basis with only 1 xi direction",err,error,*999)
5827 CALL flagerror(
"Can only set collapsed xi directions for a Lagrange Hermite tensor product basis type",err,error,*999)
5831 CALL flagerror(
"Basis is not associated",err,error,*999)
5834 exits(
"BASIS_COLLAPSED_XI_SET_PTR")
5836 999 errorsexits(
"BASIS_COLLAPSED_XI_SET_PTR",err,error)
5849 INTEGER(INTG),
INTENT(IN) :: USER_NUMBER
5851 INTEGER(INTG),
INTENT(OUT) :: ERR
5855 enters(
"BASIS_USER_NUMBER_FIND",err,error,*999)
5859 exits(
"BASIS_USER_NUMBER_FIND")
5861 999 errorsexits(
"BASIS_USER_NUMBER_FIND",err,error)
5874 INTEGER(INTG),
INTENT(IN) :: N
5875 REAL(DP),
INTENT(IN) :: ALPHA
5876 REAL(DP),
INTENT(IN) :: BETA
5877 REAL(DP),
INTENT(OUT) :: X(n)
5878 REAL(DP),
INTENT(OUT) :: W(n)
5879 INTEGER(INTG),
INTENT(OUT) :: ERR
5883 REAL(DP) :: DIFFERENCE,T1,T2
5885 INTEGER(INTG) :: GAUSS_START(4) = [ 0,1,3,6 ]
5886 REAL(DP) :: XIG(10),WIG(10)
5897 xig = [ 0.500000000000000_dp, &
5898 & (-1.0_dp/sqrt(3.0_dp)+1.0_dp)/2.0_dp,(+1.0_dp/sqrt(3.0_dp)+1.0_dp)/2.0_dp, &
5899 & (-sqrt(0.6_dp)+1.0_dp)/2.0_dp, 0.5_dp, (+sqrt(0.6_dp)+1.0_dp)/2.0_dp, &
5900 & (-sqrt((3.0_dp+2.0_dp*sqrt(6.0_dp/5.0_dp))/7.0_dp)+1.0_dp)/2.0_dp, &
5901 & (-sqrt((3.0_dp-2.0_dp*sqrt(6.0_dp/5.0_dp))/7.0_dp)+1.0_dp)/2.0_dp, &
5902 & (+sqrt((3.0_dp-2.0_dp*sqrt(6.0_dp/5.0_dp))/7.0_dp)+1.0_dp)/2.0_dp, &
5903 & (+sqrt((3.0_dp+2.0_dp*sqrt(6.0_dp/5.0_dp))/7.0_dp)+1.0_dp)/2.0_dp ]
5904 wig = [ 1.000000000000000_dp, &
5905 & 0.500000000000000_dp,0.500000000000000_dp, &
5906 & 2.5_dp/9.0_dp, 4.0_dp/9.0_dp, 2.5_dp/9.0_dp, &
5907 & (18.0_dp-sqrt(30.0_dp))/72.0_dp, &
5908 & (18.0_dp+sqrt(30.0_dp))/72.0_dp, &
5909 & (18.0_dp+sqrt(30.0_dp))/72.0_dp, &
5910 & (18.0_dp-sqrt(30.0_dp))/72.0_dp ]
5913 enters(
"GAUSS_LEGENDRE",err,error,*999)
5915 IF(n>=1.AND.n<=4)
THEN 5917 x(i)=xig(gauss_start(n)+i)
5918 w(i)=wig(gauss_start(n)+i)
5921 CALL flagerror(
"Invalid number of Gauss points. Not implemented",err,error,*999)
5926 CALL write_string_vector(
diagnostic_output_type,1,1,n,5,5,x,
'("Gauss point locations :",5(X,F13.5))',
'(23X,5(X,F13.5))', &
5928 CALL write_string_vector(
diagnostic_output_type,1,1,n,5,5,w,
'("Gauss point weights :",5(X,F13.5))',
'(23X,5(X,F13.5))', &
5935 t1=t1+((x(i)+1.0_dp)*w(i))
5937 t2=(beta**2.0_dp/2.0_dp+beta)-(alpha**2.0_dp/2.0_dp-alpha)
5938 difference=abs(t2-t1)
5944 exits(
"GAUSS_LEGENDRE")
5946 999 errorsexits(
"GAUSS_LEGENDRE",err,error)
5959 SUBROUTINE gauss_simplex(ORDER,NUMBER_OF_VERTICES,N,X,W,ERR,ERROR,*)
5962 INTEGER(INTG),
INTENT(IN) :: ORDER
5963 INTEGER(INTG),
INTENT(IN) :: NUMBER_OF_VERTICES
5964 INTEGER(INTG),
INTENT(OUT) :: N
5965 REAL(DP),
INTENT(OUT) :: X(:,:)
5966 REAL(DP),
INTENT(OUT) :: W(:)
5967 INTEGER(INTG),
INTENT(OUT) :: ERR
5971 REAL(DP) :: ALPHA_1,ALPHA_2,BETA,LAMBDA,L_C,L1_ALPHA_1,L2_ALPHA_1,L3_ALPHA_1,L4_ALPHA_1,L1_ALPHA_2,L2_ALPHA_2,L3_ALPHA_2, &
5972 & L4_ALPHA_2,L1_BETA,L2_BETA,L3_BETA,L4_BETA,W_C,W_ALPHA_1,W_ALPHA_2,W_BETA,ACOS_ARG
5975 enters(
"GAUSS_SIMPLEX",err,error,*999)
5976 IF(
SIZE(x,1)>=(number_of_vertices))
THEN 5977 SELECT CASE(number_of_vertices)
5983 IF(
SIZE(x,2)>=n)
THEN 5984 IF(
SIZE(w,1)>=n)
THEN 5985 l_c=1.0_dp/
REAL(number_of_vertices,
dp)
5992 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
5994 CALL flagerror(local_error,err,error,*999)
5997 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
5999 CALL flagerror(local_error,err,error,*999)
6003 IF(
SIZE(x,2)>=n)
THEN 6004 IF(
SIZE(w,1)>=n)
THEN 6005 alpha_1=1.0_dp/sqrt(
REAL(number_of_vertices,
dp)+1.0_DP)
6006 w_alpha_1=1.0_dp/
REAL(number_of_vertices,
dp)
6007 l1_alpha_1=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_1)/
REAL(number_of_vertices,
dp)
6008 l2_alpha_1=1.0_dp-l1_alpha_1
6018 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6020 CALL flagerror(local_error,err,error,*999)
6023 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6025 CALL flagerror(local_error,err,error,*999)
6029 IF(
SIZE(x,2)>=n)
THEN 6030 IF(
SIZE(w,1)>=n)
THEN 6031 l_c=1.0_dp/
REAL(number_of_vertices,
dp)
6032 w_c=-1.0_dp*
REAL(number_of_vertices,
dp)*
REAL(number_of_vertices,
dp)/ &
6033 & (REAL(NUMBER_OF_VERTICES,DP)*(REAL(NUMBER_OF_VERTICES,DP)+1.0_DP))
6034 ALPHA_1=2.0_dp/(
REAL(number_of_vertices,
dp)+2.0_DP)
6035 w_alpha_1=(
REAL(number_of_vertices,
dp)+2.0_DP)*(
REAL(number_of_vertices,
dp)+2.0_DP)/ &
6036 & (4.0_DP*REAL(NUMBER_OF_VERTICES,DP)*(REAL(NUMBER_OF_VERTICES,DP)+1.0_DP))
6037 L1_ALPHA_1=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_1)/
REAL(number_of_vertices,
dp)
6038 l2_alpha_1=1.0_dp-l1_alpha_1
6052 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6054 CALL flagerror(local_error,err,error,*999)
6057 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6059 CALL flagerror(local_error,err,error,*999)
6062 CALL flagerror(
"Not implemented",err,error,*999)
6063 IF(
SIZE(x,2)>=n)
THEN 6064 IF(
SIZE(w,1)>=n)
THEN 6066 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6068 CALL flagerror(local_error,err,error,*999)
6071 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6073 CALL flagerror(local_error,err,error,*999)
6076 CALL flagerror(
"Not implemented",err,error,*999)
6077 IF(
SIZE(x,2)>=n)
THEN 6078 IF(
SIZE(w,1)>=n)
THEN 6080 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6082 CALL flagerror(local_error,err,error,*999)
6085 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6087 CALL flagerror(local_error,err,error,*999)
6091 &
" is an invalid Gauss order. You must have an order between 1 and 5" 6092 CALL flagerror(local_error,err,error,*999)
6099 IF(
SIZE(x,2)>=n)
THEN 6100 IF(
SIZE(w,1)>=n)
THEN 6101 l_c=1.0_dp/
REAL(number_of_vertices,
dp)
6109 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6111 CALL flagerror(local_error,err,error,*999)
6114 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6116 CALL flagerror(local_error,err,error,*999)
6120 IF(
SIZE(x,2)>=n)
THEN 6121 IF(
SIZE(w,1)>=n)
THEN 6122 alpha_1=-1.0_dp/sqrt(
REAL(number_of_vertices,
dp)+1.0_DP)
6123 w_alpha_1=1.0_dp/
REAL(number_of_vertices,
dp)
6124 l1_alpha_1=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_1)/
REAL(number_of_vertices,
dp)
6125 l2_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6126 l3_alpha_1=1.0_dp-l1_alpha_1-l2_alpha_1
6131 w(1)=w_alpha_1/2.0_dp
6136 w(2)=w_alpha_1/2.0_dp
6141 w(3)=w_alpha_1/2.0_dp
6143 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6145 CALL flagerror(local_error,err,error,*999)
6148 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6150 CALL flagerror(local_error,err,error,*999)
6154 IF(
SIZE(x,2)>=n)
THEN 6155 IF(
SIZE(w,1)>=n)
THEN 6156 l_c=1.0_dp/
REAL(number_of_vertices,
dp)
6157 w_c=-1.0_dp*
REAL(number_of_vertices,
dp)*
REAL(number_of_vertices,
dp)/ &
6158 & (REAL(NUMBER_OF_VERTICES,DP)*(REAL(NUMBER_OF_VERTICES,DP)+1.0_DP))
6159 ALPHA_1=2.0_dp/(
REAL(number_of_vertices,
dp)+2.0_DP)
6160 w_alpha_1=(
REAL(number_of_vertices,
dp)+2.0_DP)*(
REAL(number_of_vertices,
dp)+2.0_DP)/ &
6161 & (4.0_DP*REAL(NUMBER_OF_VERTICES,DP)*(REAL(NUMBER_OF_VERTICES,DP)+1.0_DP))
6162 L1_ALPHA_1=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_1)/
REAL(number_of_vertices,
dp)
6163 l2_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6164 l3_alpha_1=1.0_dp-l1_alpha_1-l2_alpha_1
6174 w(2)=w_alpha_1/2.0_dp
6179 w(3)=w_alpha_1/2.0_dp
6184 w(4)=w_alpha_1/2.0_dp
6186 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6188 CALL flagerror(local_error,err,error,*999)
6191 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6193 CALL flagerror(local_error,err,error,*999)
6197 IF(
SIZE(x,2)>=n)
THEN 6198 IF(
SIZE(w,1)>=n)
THEN 6199 alpha_1=(-10.0_dp+5.0_dp*sqrt(10.0_dp)+sqrt(950.0_dp-220.0_dp*sqrt(10.0_dp)))/30.0_dp
6200 alpha_2=(-10.0_dp+5.0_dp*sqrt(10.0_dp)-sqrt(950.0_dp-220.0_dp*sqrt(10.0_dp)))/30.0_dp
6201 w_alpha_1=(5.0_dp*alpha_2-2.0_dp)/(60.0_dp*alpha_1*alpha_1*(alpha_2-alpha_1))
6202 w_alpha_2=(5.0_dp*alpha_1-2.0_dp)/(60.0_dp*alpha_2*alpha_2*(alpha_1-alpha_2))
6203 l1_alpha_1=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_1)/
REAL(number_of_vertices,
dp)
6204 l2_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6205 l3_alpha_1=1.0_dp-l1_alpha_1-l2_alpha_1
6206 l1_alpha_2=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_2)/
REAL(number_of_vertices,
dp)
6207 l2_alpha_2=(1.0_dp-alpha_2)/
REAL(number_of_vertices,
dp)
6208 l3_alpha_2=1.0_dp-l1_alpha_2-l2_alpha_2
6213 w(1)=w_alpha_1/2.0_dp
6218 w(2)=w_alpha_1/2.0_dp
6223 w(3)=w_alpha_1/2.0_dp
6228 w(4)=w_alpha_2/2.0_dp
6233 w(5)=w_alpha_2/2.0_dp
6238 w(6)=w_alpha_2/2.0_dp
6240 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6242 CALL flagerror(local_error,err,error,*999)
6245 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6247 CALL flagerror(local_error,err,error,*999)
6251 IF(
SIZE(x,2)>=n)
THEN 6252 IF(
SIZE(w,1)>=n)
THEN 6253 l_c=1.0_dp/
REAL(number_of_vertices,
dp)
6255 alpha_1=(1.0_dp+sqrt(15.0_dp))/7.0_dp
6256 alpha_2=(1.0_dp-sqrt(15.0_dp))/7.0_dp
6257 w_alpha_1=(155.0_dp-sqrt(15.0_dp))/1200.0_dp
6258 w_alpha_2=(155.0_dp+sqrt(15.0_dp))/1200.0_dp
6259 l1_alpha_1=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_1)/
REAL(number_of_vertices,
dp)
6260 l2_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6261 l3_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6262 l4_alpha_1=1.0_dp-l1_alpha_1-l2_alpha_1-l3_alpha_1
6263 l1_alpha_2=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_2)/
REAL(number_of_vertices,
dp)
6264 l2_alpha_2=(1.0_dp-alpha_2)/
REAL(number_of_vertices,
dp)
6265 l3_alpha_2=(1.0_dp-alpha_2)/
REAL(number_of_vertices,
dp)
6266 l4_alpha_2=1.0_dp-l1_alpha_2-l2_alpha_2-l3_alpha_2
6276 w(2)=w_alpha_1/2.0_dp
6281 w(3)=w_alpha_1/2.0_dp
6286 w(4)=w_alpha_1/2.0_dp
6291 w(5)=w_alpha_2/2.0_dp
6296 w(6)=w_alpha_2/2.0_dp
6301 w(7)=w_alpha_2/2.0_dp
6303 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6305 CALL flagerror(local_error,err,error,*999)
6308 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6310 CALL flagerror(local_error,err,error,*999)
6314 &
" is an invalid Gauss order. You must have an order between 1 and 5" 6315 CALL flagerror(local_error,err,error,*999)
6322 IF(
SIZE(x,2)>=n)
THEN 6323 IF(
SIZE(w,1)>=n)
THEN 6324 l_c=1.0_dp/
REAL(number_of_vertices,
dp)
6332 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6334 CALL flagerror(local_error,err,error,*999)
6337 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6339 CALL flagerror(local_error,err,error,*999)
6343 IF(
SIZE(x,2)>=n)
THEN 6344 IF(
SIZE(w,1)>=n)
THEN 6345 alpha_1=1.0_dp/sqrt(
REAL(number_of_vertices,
dp)+1.0_DP)
6346 w_alpha_1=1.0_dp/
REAL(number_of_vertices,
dp)
6347 l1_alpha_1=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_1)/
REAL(number_of_vertices,
dp)
6348 l2_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6349 l3_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6350 l4_alpha_1=1.0_dp-l1_alpha_1-l2_alpha_1-l3_alpha_1
6356 w(1)=w_alpha_1/6.0_dp
6362 w(2)=w_alpha_1/6.0_dp
6368 w(3)=w_alpha_1/6.0_dp
6374 w(4)=w_alpha_1/6.0_dp
6376 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6378 CALL flagerror(local_error,err,error,*999)
6381 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6383 CALL flagerror(local_error,err,error,*999)
6387 IF(
SIZE(x,2)>=n)
THEN 6388 IF(
SIZE(w,1)>=n)
THEN 6389 l_c=1.0_dp/
REAL(number_of_vertices,
dp)
6390 w_c=-1.0_dp*
REAL(number_of_vertices,
dp)*
REAL(number_of_vertices,
dp)/ &
6391 & (REAL(NUMBER_OF_VERTICES,DP)*(REAL(NUMBER_OF_VERTICES,DP)+1.0_DP))
6392 ALPHA_1=2.0_dp/(
REAL(number_of_vertices,
dp)+2.0_DP)
6393 w_alpha_1=(
REAL(number_of_vertices,
dp)+2.0_DP)*(
REAL(number_of_vertices,
dp)+2.0_DP)/ &
6394 & (4.0_DP*REAL(NUMBER_OF_VERTICES,DP)*(REAL(NUMBER_OF_VERTICES,DP)+1.0_DP))
6395 L1_ALPHA_1=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_1)/
REAL(number_of_vertices,
dp)
6396 l2_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6397 l3_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6398 l4_alpha_1=1.0_dp-l1_alpha_1-l2_alpha_1-l3_alpha_1
6410 w(2)=w_alpha_1/6.0_dp
6416 w(3)=w_alpha_1/6.0_dp
6422 w(4)=w_alpha_1/6.0_dp
6428 w(5)=w_alpha_1/6.0_dp
6430 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6432 CALL flagerror(local_error,err,error,*999)
6435 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6437 CALL flagerror(local_error,err,error,*999)
6441 IF(
SIZE(x,2)>=n)
THEN 6442 IF(
SIZE(w,1)>=n)
THEN 6443 l_c=1.0_dp/
REAL(number_of_vertices,
dp)
6444 w_c=-148.0_dp/1875.0_dp
6445 alpha_1=5.0_dp/7.0_dp
6446 beta=sqrt(70.0_dp)/28.0_dp
6447 w_alpha_1=343.0_dp/7500.0_dp
6448 w_beta=56.0_dp/375.0_dp
6449 l1_alpha_1=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_1)/
REAL(number_of_vertices,
dp)
6450 l2_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6451 l3_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6452 l4_alpha_1=1.0_dp-l1_alpha_1-l2_alpha_1-l3_alpha_1
6453 l1_beta=(1.0_dp+(
REAL(number_of_vertices,
dp)-2.0_DP)*beta)/
REAL(number_of_vertices,
dp)
6455 l3_beta=(1.0_dp-2.0_dp*beta)/
REAL(number_of_vertices,
dp)
6456 l4_beta=1.0_dp-l1_beta-l2_beta-l3_beta
6468 w(2)=w_alpha_1/6.0_dp
6474 w(3)=w_alpha_1/6.0_dp
6480 w(4)=w_alpha_1/6.0_dp
6486 w(5)=w_alpha_1/6.0_dp
6524 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6526 CALL flagerror(local_error,err,error,*999)
6529 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6531 CALL flagerror(local_error,err,error,*999)
6535 IF(
SIZE(x,2)>=n)
THEN 6536 IF(
SIZE(w,1)>=n)
THEN 6537 acos_arg=67.0_dp*sqrt(79.0_dp)/24964.0_dp
6539 lambda=4.0_dp/27.0_dp*(4.0_dp*sqrt(79.0_dp)*cos(((acos(acos_arg)+
twopi)/3.0_dp))+71.0_dp)
6541 alpha_1=(sqrt(9.0_dp*lambda*lambda-248.0_dp*lambda+1680.0_dp)+28.0_dp-3.0_dp*lambda)/ &
6542 & (112.0_dp-10.0_dp*lambda)
6543 alpha_2=(-1.0_dp*sqrt(9.0_dp*lambda*lambda-248.0_dp*lambda+1680.0_dp)+28.0_dp-3.0_dp*lambda)/ &
6544 & (112.0_dp-10.0_dp*lambda)
6545 beta=1.0_dp/sqrt(lambda)
6546 w_alpha_1=((21.0_dp-lambda)*alpha_2-7.0_dp)/(420.0_dp*alpha_1*alpha_1*(alpha_2-alpha_1))
6547 w_alpha_2=((21.0_dp-lambda)*alpha_1-7.0_dp)/(420.0_dp*alpha_2*alpha_2*(alpha_1-alpha_2))
6548 w_beta=lambda*lambda/840.0_dp
6549 l1_alpha_1=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_1)/
REAL(number_of_vertices,
dp)
6550 l2_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6551 l3_alpha_1=(1.0_dp-alpha_1)/
REAL(number_of_vertices,
dp)
6552 l4_alpha_1=1.0_dp-l1_alpha_1-l2_alpha_1-l3_alpha_1
6553 l1_alpha_2=(1.0_dp+(
REAL(number_of_vertices,
dp)-1.0_DP)*alpha_2)/
REAL(number_of_vertices,
dp)
6554 l2_alpha_2=(1.0_dp-alpha_2)/
REAL(number_of_vertices,
dp)
6555 l3_alpha_2=(1.0_dp-alpha_2)/
REAL(number_of_vertices,
dp)
6556 l4_alpha_2=1.0_dp-l1_alpha_2-l2_alpha_2-l3_alpha_2
6557 l1_beta=(1.0_dp+(
REAL(number_of_vertices,
dp)-2.0_DP)*beta)/
REAL(number_of_vertices,
dp)
6558 l2_beta=(1.0_dp+(
REAL(number_of_vertices,
dp)-2.0_DP)*beta)/
REAL(number_of_vertices,
dp)
6559 l3_beta=(1.0_dp-2.0_dp*beta)/
REAL(number_of_vertices,
dp)
6560 l4_beta=1.0_dp-l1_beta-l2_beta-l3_beta
6566 w(1)=w_alpha_1/6.0_dp
6572 w(2)=w_alpha_1/6.0_dp
6578 w(3)=w_alpha_1/6.0_dp
6584 w(4)=w_alpha_1/6.0_dp
6590 w(5)=w_alpha_2/6.0_dp
6596 w(6)=w_alpha_2/6.0_dp
6602 w(7)=w_alpha_2/6.0_dp
6608 w(8)=w_alpha_2/6.0_dp
6646 local_error=
"The first dimension of the W array is "//
trim(
number_to_vstring(
SIZE(w,1),
"*",err,error))// &
6648 CALL flagerror(local_error,err,error,*999)
6651 local_error=
"The second dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,2),
"*",err,error))// &
6653 CALL flagerror(local_error,err,error,*999)
6657 &
" is an invalid Gauss order. You must have an order between 1 and 5" 6658 CALL flagerror(local_error,err,error,*999)
6662 &
" is an invalid number of vertices. You must have between 2 and 4 vertices" 6663 CALL flagerror(local_error,err,error,*999)
6666 local_error=
"The first dimension of the X array is "//
trim(
number_to_vstring(
SIZE(x,1),
"*",err,error))// &
6667 &
" and it must be >= the number of vertices" 6668 CALL flagerror(local_error,err,error,*999)
6678 CALL write_string_vector(
diagnostic_output_type,1,1,order,4,4,x(:,ng),
'(" Location(nic) :",4(X,F13.5))', &
6679 &
'(23X,4(X,F13.5))',err,error,*999)
6687 exits(
"GAUSS_SIMPLEX")
6689 999 errorsexits(
"GAUSS_SIMPLEX",err,error)
6701 INTEGER(INTG),
INTENT(IN) :: NODE_INDEX
6702 INTEGER(INTG),
INTENT(IN) :: NODE_DERIVATIVE_INDEX
6703 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIVATIVE_INDEX
6704 REAL(DP),
INTENT(IN) :: XI
6705 INTEGER(INTG),
INTENT(OUT) :: ERR
6708 REAL(DP) :: HERMITE_CUBIC_EVALUATE
6711 enters(
"HERMITE_CUBIC_EVALUATE",err,error,*999)
6713 hermite_cubic_evaluate=0.0_dp
6714 SELECT CASE(partial_derivative_index)
6716 SELECT CASE(node_index)
6718 SELECT CASE(node_derivative_index)
6720 hermite_cubic_evaluate=(2.0_dp*xi-3.0_dp)*xi*xi+1.0_dp
6722 hermite_cubic_evaluate=((xi-2.0_dp)*xi+1.0_dp)*xi
6724 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6727 SELECT CASE(node_derivative_index)
6729 hermite_cubic_evaluate=xi*xi*(3.0_dp-2.0_dp*xi)
6731 hermite_cubic_evaluate=xi*xi*(xi-1.0_dp)
6733 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6736 CALL flagerror(
"Invalid node index",err,error,*999)
6739 SELECT CASE(node_index)
6741 SELECT CASE(node_derivative_index)
6743 hermite_cubic_evaluate=6.0_dp*xi*(xi-1.0_dp)
6745 hermite_cubic_evaluate=(3.0_dp*xi-4.0_dp)*xi+1.0_dp
6747 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6750 SELECT CASE(node_derivative_index)
6752 hermite_cubic_evaluate=6.0_dp*xi*(1.0_dp-xi)
6754 hermite_cubic_evaluate=xi*(3.0_dp*xi-2.0_dp)
6756 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6759 CALL flagerror(
"Invalid node index",err,error,*999)
6762 SELECT CASE(node_index)
6764 SELECT CASE(node_derivative_index)
6766 hermite_cubic_evaluate=12.0_dp*xi-6.0_dp
6768 hermite_cubic_evaluate=6.0_dp*xi-4.0_dp
6770 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6773 SELECT CASE(node_derivative_index)
6775 hermite_cubic_evaluate=6.0_dp-12.0_dp*xi
6777 hermite_cubic_evaluate=6.0_dp*xi-2.0_dp
6779 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6782 CALL flagerror(
"Invalid node index",err,error,*999)
6785 CALL flagerror(
"Invalid partial derivative index",err,error,*999)
6788 exits(
"HERMITE_CUBIC_EVALUATE")
6790 999 errorsexits(
"HERMITE_CUBIC_EVALUATE",err,error)
6809 FUNCTION hermite_quadratic_evaluate(NODE_INDEX,NODE_DERIVATIVE_INDEX,PARTIAL_DERIVATIVE_INDEX,SPECIAL_NODE_INDEX,XI,ERR,ERROR)
6812 INTEGER(INTG),
INTENT(IN) :: NODE_INDEX
6813 INTEGER(INTG),
INTENT(IN) :: NODE_DERIVATIVE_INDEX
6814 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIVATIVE_INDEX
6815 INTEGER(INTG),
INTENT(IN) :: SPECIAL_NODE_INDEX
6816 REAL(DP),
INTENT(IN) :: XI
6817 INTEGER(INTG),
INTENT(OUT) :: ERR
6820 REAL(DP) :: HERMITE_QUADRATIC_EVALUATE
6823 enters(
"HERMITE_QUADRATIC_EVALUATE",err,error,*999)
6825 hermite_quadratic_evaluate=0.0_dp
6826 SELECT CASE(special_node_index)
6828 SELECT CASE(partial_derivative_index)
6830 SELECT CASE(node_index)
6832 SELECT CASE(node_derivative_index)
6834 hermite_quadratic_evaluate=(xi-2.0_dp)*xi+1.0_dp
6836 hermite_quadratic_evaluate=0.0_dp
6838 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6841 SELECT CASE(node_derivative_index)
6843 hermite_quadratic_evaluate=(2.0_dp-xi)*xi
6845 hermite_quadratic_evaluate=(xi-1.0_dp)*xi
6847 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6850 CALL flagerror(
"Invalid node index",err,error,*999)
6853 SELECT CASE(node_index)
6855 SELECT CASE(node_derivative_index)
6857 hermite_quadratic_evaluate=2.0_dp*xi-2.0_dp
6859 hermite_quadratic_evaluate=0.0_dp
6861 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6864 SELECT CASE(node_derivative_index)
6866 hermite_quadratic_evaluate=-2.0_dp*xi+2.0_dp
6868 hermite_quadratic_evaluate=2.0_dp*xi-1.0_dp
6870 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6873 CALL flagerror(
"Invalid node index",err,error,*999)
6876 SELECT CASE(node_index)
6878 SELECT CASE(node_derivative_index)
6880 hermite_quadratic_evaluate=2.0_dp
6882 hermite_quadratic_evaluate=0.0_dp
6884 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6887 SELECT CASE(node_derivative_index)
6889 hermite_quadratic_evaluate=-2.0_dp
6891 hermite_quadratic_evaluate=2.0_dp
6893 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6896 CALL flagerror(
"Invalid node index",err,error,*999)
6899 CALL flagerror(
"Invalid partial derivative index",err,error,*999)
6902 SELECT CASE(partial_derivative_index)
6904 SELECT CASE(node_index)
6906 SELECT CASE(node_derivative_index)
6908 hermite_quadratic_evaluate=1.0_dp-xi*xi
6910 hermite_quadratic_evaluate=xi*(1.0_dp-xi)
6912 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6915 SELECT CASE(node_derivative_index)
6917 hermite_quadratic_evaluate=xi*xi
6919 hermite_quadratic_evaluate=0.0_dp
6921 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6924 CALL flagerror(
"Invalid node index",err,error,*999)
6927 SELECT CASE(node_index)
6929 SELECT CASE(node_derivative_index)
6931 hermite_quadratic_evaluate=-2.0_dp*xi
6933 hermite_quadratic_evaluate=1.0_dp-2.0_dp*xi
6935 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6938 SELECT CASE(node_derivative_index)
6940 hermite_quadratic_evaluate=2.0_dp*xi
6942 hermite_quadratic_evaluate=0.0_dp
6944 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6947 CALL flagerror(
"Invalid node index",err,error,*999)
6950 SELECT CASE(node_index)
6952 SELECT CASE(node_derivative_index)
6954 hermite_quadratic_evaluate=-2.0_dp
6956 hermite_quadratic_evaluate=-2.0_dp
6958 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6961 SELECT CASE(node_derivative_index)
6963 hermite_quadratic_evaluate=2.0_dp
6965 hermite_quadratic_evaluate=0.0_dp
6967 CALL flagerror(
"Invalid node derivative index",err,error,*999)
6970 CALL flagerror(
"Invalid node index",err,error,*999)
6973 CALL flagerror(
"Invalid partial derivative index",err,error,*999)
6976 CALL flagerror(
"Invalid special node index",err,error,*999)
6979 exits(
"HERMITE_QUADRATIC_EVALUATE")
6981 999 errorsexits(
"HERMITE_QUADRATIC_EVALUATE",err,error)
6993 INTEGER(INTG),
INTENT(IN) :: NODE_INDEX
6994 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIVATIVE_INDEX
6995 REAL(DP),
INTENT(IN) :: XI
6996 INTEGER(INTG),
INTENT(OUT) :: ERR
6999 REAL(DP) :: LAGRANGE_CUBIC_EVALUATE
7002 enters(
"LAGRANGE_CUBIC_EVALUATE",err,error,*999)
7004 lagrange_cubic_evaluate=0.0_dp
7005 SELECT CASE(partial_derivative_index)
7007 SELECT CASE(node_index)
7009 lagrange_cubic_evaluate=0.5_dp*(3.0_dp*xi-1.0_dp)*(3.0_dp*xi-2.0_dp)*(1.0_dp-xi)
7011 lagrange_cubic_evaluate=4.5_dp*xi*(3.0_dp*xi-2.0_dp)*(xi-1.0_dp)
7013 lagrange_cubic_evaluate=4.5_dp*xi*(3.0_dp*xi-1.0_dp)*(1.0_dp-xi)
7015 lagrange_cubic_evaluate=0.5_dp*xi*(3.0_dp*xi-1.0_dp)*(3.0_dp*xi-2.0_dp)
7017 CALL flagerror(
"Invalid node index",err,error,*999)
7020 SELECT CASE(node_index)
7022 lagrange_cubic_evaluate=-13.5_dp*xi*xi+18.0_dp*xi-5.5_dp
7024 lagrange_cubic_evaluate= 40.5_dp*xi*xi-45.0_dp*xi+9.0_dp
7026 lagrange_cubic_evaluate=-40.5_dp*xi*xi+36.0_dp*xi-4.5_dp
7028 lagrange_cubic_evaluate= 13.5_dp*xi*xi- 9.0_dp*xi+1.0_dp
7030 CALL flagerror(
"Invalid node index",err,error,*999)
7033 SELECT CASE(node_index)
7035 lagrange_cubic_evaluate=9.0_dp*(2.0_dp-3.0_dp*xi)
7037 lagrange_cubic_evaluate=9.0_dp*(9.0_dp*xi-5.0_dp)
7039 lagrange_cubic_evaluate=9.0_dp*(4.0_dp-9.0_dp*xi)
7041 lagrange_cubic_evaluate=9.0_dp*(3.0_dp*xi-1.0_dp)
7043 CALL flagerror(
"Invalid node index",err,error,*999)
7046 CALL flagerror(
"Invalid partial derivative index",err,error,*999)
7049 exits(
"LAGRANGE_CUBIC_EVALUATE")
7051 999 errorsexits(
"LAGRANGE_CUBIC_EVALUATE",err,error)
7063 INTEGER(INTG),
INTENT(IN) :: NODE_INDEX
7064 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIVATIVE_INDEX
7065 REAL(DP),
INTENT(IN) :: XI
7066 INTEGER(INTG),
INTENT(OUT) :: ERR
7069 REAL(DP) :: LAGRANGE_LINEAR_EVALUATE
7072 enters(
"LAGRANGE_LINEAR_EVALUATE",err,error,*999)
7074 lagrange_linear_evaluate=0.0_dp
7075 SELECT CASE(partial_derivative_index)
7077 SELECT CASE(node_index)
7079 lagrange_linear_evaluate=1.0_dp-xi
7081 lagrange_linear_evaluate=xi
7083 CALL flagerror(
"Invalid node index",err,error,*999)
7086 SELECT CASE(node_index)
7088 lagrange_linear_evaluate=-1.0_dp
7090 lagrange_linear_evaluate=1.0_dp
7092 CALL flagerror(
"Invalid node index",err,error,*999)
7095 SELECT CASE(node_index)
7097 lagrange_linear_evaluate=0.0_dp
7099 lagrange_linear_evaluate=0.0_dp
7101 CALL flagerror(
"Invalid node index",err,error,*999)
7104 CALL flagerror(
"Invalid partial derivative index",err,error,*999)
7107 exits(
"LAGRANGE_LINEAR_EVALUATE")
7109 999 errorsexits(
"LAGRANGE_LINEAR_EVALUATE",err,error)
7121 INTEGER(INTG),
INTENT(IN) :: NODE_INDEX
7122 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIVATIVE_INDEX
7123 REAL(DP),
INTENT(IN) :: XI
7124 INTEGER(INTG),
INTENT(OUT) :: ERR
7127 REAL(DP) :: LAGRANGE_QUADRATIC_EVALUATE
7130 enters(
"LAGRANGE_QUADRATIC_EVALUATE",err,error,*999)
7132 lagrange_quadratic_evaluate=0.0_dp
7133 SELECT CASE(partial_derivative_index)
7135 SELECT CASE(node_index)
7137 lagrange_quadratic_evaluate=1.0_dp-3.0_dp*xi+2.0_dp*xi*xi
7139 lagrange_quadratic_evaluate=4.0_dp*xi*(1.0_dp-xi)
7141 lagrange_quadratic_evaluate=xi*(xi+xi-1.0_dp)
7143 CALL flagerror(
"Invalid node index",err,error,*999)
7146 SELECT CASE(node_index)
7148 lagrange_quadratic_evaluate=4.0_dp*xi-3.0_dp
7150 lagrange_quadratic_evaluate=4.0_dp-8.0_dp*xi
7152 lagrange_quadratic_evaluate=4.0_dp*xi-1.0_dp
7154 CALL flagerror(
"Invalid node index",err,error,*999)
7157 SELECT CASE(node_index)
7159 lagrange_quadratic_evaluate=4.0_dp
7161 lagrange_quadratic_evaluate=-8.0_dp
7163 lagrange_quadratic_evaluate=4.0_dp
7165 CALL flagerror(
"Invalid node index",err,error,*999)
7168 CALL flagerror(
"Invalid partial derivative index",err,error,*999)
7171 exits(
"LAGRANGE_QUADRATIC_EVALUATE")
7173 999 errorsexits(
"LAGRANGE_QUADRATIC_EVALUATE",err,error)
7186 INTEGER(INTG),
INTENT(IN) :: NODE_INDEX
7187 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIVATIVE_INDEX
7188 REAL(DP),
INTENT(IN) :: XL
7189 INTEGER(INTG),
INTENT(OUT) :: ERR
7192 REAL(DP) :: SIMPLEX_CUBIC_EVALUATE_DP
7195 enters(
"SIMPLEX_CUBIC_EVALUATE_DP",err,error,*999)
7197 simplex_cubic_evaluate_dp=0.0_dp
7199 SELECT CASE(partial_derivative_index)
7201 SELECT CASE(node_index)
7203 simplex_cubic_evaluate_dp=1.0
7205 simplex_cubic_evaluate_dp=3.0_dp*xl
7207 simplex_cubic_evaluate_dp=3.0_dp/2.0_dp*xl*(3.0_dp*xl-1.0_dp)
7209 simplex_cubic_evaluate_dp=0.5_dp*xl*(3.0_dp*xl-1.0_dp)*(3.0_dp*xl-2.0_dp)
7211 CALL flagerror(
"Invalid node index.",err,error,*999)
7214 SELECT CASE(node_index)
7216 simplex_cubic_evaluate_dp=0.0_dp
7218 simplex_cubic_evaluate_dp=3.0_dp
7220 simplex_cubic_evaluate_dp=3.0_dp/2.0_dp*(6.0_dp*xl-1)
7222 simplex_cubic_evaluate_dp=13.5_dp*xl*xl-9.0_dp*xl+1.0_dp
7224 CALL flagerror(
"Invalid node index.",err,error,*999)
7227 SELECT CASE(node_index)
7229 simplex_cubic_evaluate_dp=0.0_dp
7231 simplex_cubic_evaluate_dp=0.0_dp
7233 simplex_cubic_evaluate_dp=9.0_dp
7235 simplex_cubic_evaluate_dp=2.0_dp*xl-9.0_dp
7237 CALL flagerror(
"Invalid node index.",err,error,*999)
7240 SELECT CASE(node_index)
7242 simplex_cubic_evaluate_dp=0.0_dp
7244 simplex_cubic_evaluate_dp=0.0_dp
7246 simplex_cubic_evaluate_dp=0.0_dp
7248 simplex_cubic_evaluate_dp=2.0_dp
7250 CALL flagerror(
"Invalid node index.",err,error,*999)
7253 CALL flagerror(
"Invalid partial derivative index.",err,error,*999)
7256 exits(
"SIMPLEX_CUBIC_EVALUATE_DP")
7258 999 errorsexits(
"SIMPLEX_CUBIC_EVALUATE_DP",err,error)
7271 INTEGER(INTG),
INTENT(IN) :: NODE_INDEX
7272 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIVATIVE_INDEX
7273 REAL(DP),
INTENT(IN) :: XL
7274 INTEGER(INTG),
INTENT(OUT) :: ERR
7277 REAL(DP) :: SIMPLEX_LINEAR_EVALUATE_DP
7280 enters(
"SIMPLEX_LINEAR_EVALUATE_DP",err,error,*999)
7282 simplex_linear_evaluate_dp=0.0_dp
7283 SELECT CASE(partial_derivative_index)
7285 SELECT CASE(node_index)
7287 simplex_linear_evaluate_dp=1.0
7289 simplex_linear_evaluate_dp=xl
7291 CALL flagerror(
"Invalid node index",err,error,*999)
7294 SELECT CASE(node_index)
7296 simplex_linear_evaluate_dp=0.0_dp
7298 simplex_linear_evaluate_dp=1.0_dp
7300 CALL flagerror(
"Invalid node index",err,error,*999)
7303 SELECT CASE(node_index)
7305 simplex_linear_evaluate_dp=0.0_dp
7307 simplex_linear_evaluate_dp=0.0_dp
7309 CALL flagerror(
"Invalid node index",err,error,*999)
7312 SELECT CASE(node_index)
7314 simplex_linear_evaluate_dp=0.0_dp
7316 simplex_linear_evaluate_dp=0.0_dp
7318 CALL flagerror(
"Invalid node index",err,error,*999)
7321 CALL flagerror(
"Invalid partial derivative index",err,error,*999)
7324 exits(
"SIMPLEX_LINEAR_EVALUATE_DP")
7326 999 errorsexits(
"SIMPLEX_LINEAR_EVALUATE_DP",err,error)
7339 INTEGER(INTG),
INTENT(IN) :: NODE_INDEX
7340 INTEGER(INTG),
INTENT(IN) :: PARTIAL_DERIVATIVE_INDEX
7341 REAL(DP),
INTENT(IN) :: XL
7342 INTEGER(INTG),
INTENT(OUT) :: ERR
7345 REAL(DP) :: SIMPLEX_QUADRATIC_EVALUATE_DP
7348 enters(
"SIMPLEX_QUADRATIC_EVALUATE_DP",err,error,*999)
7350 simplex_quadratic_evaluate_dp=0.0_dp
7351 SELECT CASE(partial_derivative_index)
7353 SELECT CASE(node_index)
7355 simplex_quadratic_evaluate_dp=1.0_dp
7357 simplex_quadratic_evaluate_dp=2.0_dp*xl
7359 simplex_quadratic_evaluate_dp=xl*(2.0_dp*xl-1.0_dp)
7361 CALL flagerror(
"Invalid node index.",err,error,*999)
7364 SELECT CASE(node_index)
7366 simplex_quadratic_evaluate_dp=0.0_dp
7368 simplex_quadratic_evaluate_dp=2.0_dp
7370 simplex_quadratic_evaluate_dp=4.0_dp*xl-1.0_dp
7372 CALL flagerror(
"Invalid node index",err,error,*999)
7375 SELECT CASE(node_index)
7377 simplex_quadratic_evaluate_dp=0.0_dp
7379 simplex_quadratic_evaluate_dp=0.0_dp
7381 simplex_quadratic_evaluate_dp=4.0_dp
7383 CALL flagerror(
"Invalid node index.",err,error,*999)
7386 SELECT CASE(node_index)
7388 simplex_quadratic_evaluate_dp=0.0_dp
7390 simplex_quadratic_evaluate_dp=0.0_dp
7392 simplex_quadratic_evaluate_dp=0.0_dp
7394 CALL flagerror(
"Invalid node index.",err,error,*999)
7397 CALL flagerror(
"Invalid partial derivative index.",err,error,*999)
7400 exits(
"SIMPLEX_QUADRATIC_EVALUATE_DP")
7402 999 errorsexits(
"SIMPLEX_QUADRATIC_EVALUATE_DP",err,error)
This module contains all basis function routines.
Sets/changes the number of Xi directions for a basis.
real(dp) function basis_interpolate_local_face_gauss_dp(BASIS, PARTIAL_DERIV_INDEX, QUADRATURE_SCHEME, LOCAL_FACE_NUMBER, GAUSS_POINT_NUMBER, FACE_PARAMETERS, ERR, ERROR)
Interpolates the appropriate partial derivative index of the element local face parameters at a face ...
subroutine, public enters(NAME, ERR, ERROR,)
Records the entry into the named procedure and initialises the error code.
integer(intg), parameter, public basis_xi_collapsed
The Xi direction is collapsed.
integer(intg), parameter second_part_deriv
Second partial derivative i.e., d^2u/ds^2.
subroutine, public basis_user_number_find(USER_NUMBER, BASIS, ERR, ERROR,)
Finds and returns in BASIS a pointer to the basis with the number given in USER_NUMBER. If no basis with that number exits BASIS is left nullified.
integer(intg), parameter, public basis_quadratic_lagrange_interpolation
Quadratic Lagrange interpolation specification.
integer(intg), dimension(2) other_xi_directions2
OTHER_XI_DIRECTIONS2(ni) gives the other xi direction for direction ni for a two dimensional element...
subroutine, public basis_create_finish(BASIS, ERR, ERROR,)
Finishes the creation of a new basis.
subroutine basis_lhtpbasiscreate(basis, err, error,)
Creates and initialises a Lagrange-Hermite tensor product basis that has already been allocated BASIS...
Converts a number to its equivalent varying string representation.
subroutine basis_number_of_xi_set_ptr(BASIS, NUMBER_OF_XI, ERR, ERROR,)
Sets/changes the number of xi directions for a basis identified by a pointer.
integer(intg), parameter, public basis_quadratic2_hermite_interpolation
Quadratic Hermite (no derivative at xi=1) interpolation specification.
Evaluates the Lagrange/Hermite/Fourier tensor product basis function for the given basis...
Evaluates the appropriate partial derivative index for the specificied basis function at a Xi locatio...
subroutine basis_simplex_family_create(BASIS, ERR, ERROR,)
Creates and initialises a simplex basis family that has already been allocated by BASIS_CREATE_START...
subroutine, public basis_quadrature_multiple_gauss_xi_get(BASIS, SCHEME, GAUSS_POINTS, GAUSS_XI, ERR, ERROR,)
Returns the xi positions of Gauss points on a basis quadrature identified by a pointer. If no Gauss points are specified then xi positions of all Gauss points are returned.
Sets/changes the interpolation type in each Xi direction for a basis.
subroutine basis_radial_family_create(BASIS, ERR, ERROR,)
Creates and initialises a Radial basis family that has already been allocated by BASIS_CREATE_START.
subroutine basis_interpolation_xi_set_number(USER_NUMBER, INTERPOLATION_XI, ERR, ERROR,)
Sets/changes the interpolation type in each xi directions where the basis is identified by user numbe...
integer(intg), parameter no_part_deriv
No partial derivative i.e., u.
integer(intg), parameter basis_number_of_quadrature_scheme_types
The number of currently defined quadrature schemes.
integer(intg), parameter, public basis_collapsed_at_xi0
The Xi direction at the xi=0 end of this Xi direction is collapsed.
real(dp) function basis_interpolate_gauss_dp(BASIS, PARTIAL_DERIV_INDEX, QUADRATURE_SCHEME, GAUSS_POINT_NUMBER, ELEMENT_PARAMETERS, ERR, ERROR)
Interpolates the appropriate partial derivative index of the element parameters at a gauss point for ...
subroutine basis_type_set_number(USER_NUMBER, TYPE, ERR, ERROR,)
Sets/changes the type for a basis is identified by a user number.
integer(intg), parameter part_deriv_s4_s4_s4
Third partial derivative in the s4 direction i.e., d^3u/ds4^3.
integer(intg), parameter basis_transition_interpolation
Transition interpolation.
integer(intg), parameter part_deriv_s4
First partial derivative in the s4 direction i.e., du/ds4.
integer(intg), parameter, public basis_adaptive_gauss_legendre_quadrature
Adaptive Gauss-Legendre quadrature.
real(dp) function basis_evaluate_xi_dp(BASIS, ELEMENT_PARAMETER_INDEX, PARTIAL_DERIV_INDEX, XI, ERR, ERROR)
Evaluates the appropriate partial derivative index at position XI for the basis for double precision ...
integer(intg), parameter basis_serendipity_interpolation
Serendipity interpolation.
subroutine, public basis_quadrature_number_of_gauss_xi_set(BASIS, NUMBER_OF_GAUSS_XI, ERR, ERROR,)
Sets/changes the number of Gauss points in each xi direction on a basis quadrature identified by a po...
Sets/changes the collapsed Xi flags for a basis.
real(dp) function simplex_quadratic_evaluate_dp(NODE_INDEX, PARTIAL_DERIVATIVE_INDEX, XL, ERR, ERROR)
Evaluates a quadratic simpelx basis function at a specificed area position and node index and with a ...
This module contains all string manipulation and transformation routines.
integer(intg), parameter basis_singular_interpolation
Singular interpolation.
integer(intg), parameter, public basis_quadratic_simplex_interpolation
Quadratic Simplex interpolation specification.
integer(intg), parameter part_deriv_s3_s4
Cross derivative in the s3 and s4 direction i.e., d^2u/ds3ds4.
integer(intg), parameter, public basis_mid_quadrature_scheme
Identifier for a mid order quadrature scheme.
subroutine, public basis_interpolation_xi_get(BASIS, INTERPOLATION_XI, ERR, ERROR,)
Gets/changes the interpolation type in each xi directions for a basis identified by a pointer...
integer(intg), parameter basis_radial_interpolation
Radial interpolation.
integer(intg), parameter first_part_deriv
First partial derivative i.e., du/ds.
real(dp) function basis_simplex_basis_derivative_evaluate(BASIS, NODE_NUMBER, PARTIAL_DERIV_INDEX, XL, ERR, ERROR)
Evaluates partial derivatives of a simplex basis function with respect to area coordinates.
subroutine, public basis_collapsed_xi_get(BASIS, COLLAPSED_XI, ERR, ERROR,)
Gets the collapsed xi flags for a basis is identified by a a pointer.
integer(intg), parameter, public basis_b_spline_tp_type
B-spline basis type.
integer(intg), parameter, public basis_simplex_type
Simplex basis type.
integer(intg), parameter, public basis_gauss_legendre_quadrature
Gauss-Legendre quadrature.
integer(intg), parameter part_deriv_s2
First partial derivative in the s2 direction i.e., du/ds2.
integer(intg), parameter part_deriv_s3_s4_s4
Cross derivative in the s3, s4 and s4 direction i.e., d^3u/ds3ds4^2.
integer(intg), parameter basis_lagrange_interpolation
Lagrange interpolation.
subroutine, public basis_quadrature_single_gauss_xi_get(BASIS, SCHEME, GAUSS_POINT, GAUSS_XI, ERR, ERROR,)
Returns the xi positions of a Gauss point on a basis quadrature identified by a pointer.
subroutine basis_type_set_ptr(BASIS, TYPE, ERR, ERROR,)
Sets/changes the type for a basis is identified by a a pointer.
integer(intg), parameter part_deriv_s4_s4
Second partial derivative in the s4 direction i.e., d^2u/ds4ds4.
integer(intg), parameter, public basis_radial_type
Radial basis typee.
integer(intg), parameter, public basis_high_quadrature_scheme
Identifier for a high order quadrature scheme.
integer(intg), parameter part_deriv_s1_s4_s4
Cross derivative in the s2, s4 and s4 direction i.e., d^3u/ds1ds4^2.
logical, save, public diagnostics2
.TRUE. if level 2 diagnostic output is active in the current routine
subroutine basis_finalise(BASIS, ERR, ERROR,)
Finalises a basis and deallocates all memory.
integer(intg), parameter, public basis_quadratic1_interpolation_order
Quadratic (no derivative at xi=0) interpolation order.
subroutine, public basis_number_of_local_nodes_get(BASIS, NUMBER_OF_LOCAL_NODES, ERR, ERROR,)
Returns the number of local nodes in the specified basis.
This module provides an iso_varying_string module, conformant to the API specified in ISO/IEC 1539-2:...
integer(intg), parameter basis_simplex_interpolation
Simplex interpolation.
subroutine basis_initialise(BASIS, ERR, ERROR,)
Initialises a basis.
integer(intg), parameter, public basis_default_quadrature_scheme
Identifier for the default quadrature scheme.
integer(intg), parameter, public basis_gaussian_radial_interpolation
Gaussian Radial interpolation specification.
This module contains all program wide constants.
integer(intg), parameter part_deriv_s1
First partial derivative in the s1 direction i.e., du/ds1.
Flags a warning to the user.
subroutine, public bases_initialise(ERR, ERROR,)
Initialises the bases.
subroutine, public basis_quadraturelocalfacegaussevaluateset(BASIS, FACE_GAUSS_EVALUATE, ERR, ERROR,)
Sets/changes the local face Gauss evaluation flag on a basis.
integer(intg), parameter, public basis_serendipity_type
Serendipity basis type.
Contains information on the quadrature to be used for integrating a basis.
real(dp) function simplex_linear_evaluate_dp(NODE_INDEX, PARTIAL_DERIVATIVE_INDEX, XL, ERR, ERROR)
Evaluates a linear simpelx basis function at a specificed area position and node index and with a giv...
integer(intg), parameter, public basis_cubic_simplex_interpolation
Cubic Simplex interpolation specification.
subroutine, public basis_quadrature_destroy(QUADRATURE, ERR, ERROR,)
Destroys a quadrature on a given basis and deallocates all memory.
subroutine basis_collapsed_xi_set_number(USER_NUMBER, COLLAPSED_XI, ERR, ERROR,)
Sets/changes the collapsed xi flags for a basis is identified by a user number.
integer(intg), parameter, public basis_linear_simplex_interpolation
Linear Simplex interpolation specification.
recursive subroutine basis_family_destroy(USER_NUMBER, FAMILY_NUMBER, ERR, ERROR,)
Destroys a basis identified by its basis user number and family number. Called from the library visib...
logical, save, public diagnostics3
.TRUE. if level 3 diagnostic output is active in the current routine
integer, parameter dp
Double precision real kind.
integer(intg), parameter, public basis_gauss_laguerre_quadrature
Gauss-Laguerre quadrature.
subroutine basis_quadrature_finalise(BASIS, ERR, ERROR,)
Finalises a quadrature on a given basis and deallocates all memory.
subroutine gauss_simplex(ORDER, NUMBER_OF_VERTICES, N, X, W, ERR, ERROR,)
This routine calculates the weights and abscissae for a Gauss quadrature scheme for simplex elements...
subroutine, public exits(NAME)
Records the exit out of the named procedure.
subroutine basis_quadrature_type_set_number(USER_NUMBER, TYPE, ERR, ERROR,)
Sets/changes the quadrature type for a basis quadrature identified by a user number.
integer(intg), parameter, public basis_extended_lagrange_tp_type
Extendend Lagrange tensor product basis type.
This module contains all type definitions in order to avoid cyclic module references.
recursive subroutine, public basis_destroy(BASIS, ERR, ERROR,)
Destroys a basis.
integer(intg), parameter, public basis_guass_hermite_quadrature
Gauss-Hermite quadrature.
subroutine basis_gauss_points_calculate_dp(basis, order, numCoords, numberGaussPoints, gaussPoints, gaussWeights, err, error,)
Calculates the gauss points and weights for a basis function of a particular order.
integer(intg), parameter part_deriv_s1_s3_s4
Cross derivative in the s1, s3 and s4 direction i.e., d^3u/ds1ds3ds4.
integer(intg), parameter part_deriv_s3
First partial derivative in the s3 direction i.e., du/ds3.
Interpolates the appropriate partial derivative index of the elements parameters for basis function a...
integer(intg), parameter, public basis_fourier_lagrange_hermite_tp_type
Fourier-Lagrange tensor product basis type.
type(basis_functions_type), public basis_functions
The tree of defined basis functions.
integer(intg), dimension(3, 3, 2) other_xi_directions3
OTHER_XI_DIRECTIONS3(ni,nii,type) gives the other xi directions for direction ni for a three dimensio...
real(dp) function basis_simplex_basis_evaluate(BASIS, NODE_NUMBER, PARTIAL_DERIV_INDEX, XL, ERR, ERROR)
Evaluates a simplex basis function and its derivatives with respect to external coordinates. For Simplex line elements there are two area coordinates which are a function of : and .The derivatives wrt to external coordinates are then given by and . For Simplex triangle elements there are three area coordinates which are a function of and : , and . The derivatives wrt to external coordinates are then given by , , , and . For Simplex tetrahedral elements there are four area coordinates which are a function of , and : , , and . The derivatives wrt to external coordinates are then given by , , , , , , , and .
This module contains all the low-level base routines e.g., all debug, control, and low-level communic...
real(dp) function simplex_cubic_evaluate_dp(NODE_INDEX, PARTIAL_DERIVATIVE_INDEX, XL, ERR, ERROR)
Evaluates a cubic simpelx basis function at a specificed area position and node index and with a give...
Interpolates the appropriate partial derivative index of the elements parameters for basis function a...
recursive subroutine basis_destroy_number(USER_NUMBER, ERR, ERROR,)
Destroys a basis identified by its basis user number.
integer(intg), parameter part_deriv_s1_s1
Second partial derivative in the s1 direction i.e., d^2u/ds1ds1.
integer(intg), parameter, public basis_not_collapsed
The Xi direction is not collapsed.
integer(intg), parameter, public basis_cubic_interpolation_order
Cubic interpolation order.
integer(intg), parameter, public basis_quadratic_interpolation_order
Quadratic interpolation order.
integer(intg), parameter part_deriv_s2_s3
Cross derivative in the s2 and s3 direction i.e., d^2u/ds2ds3.
integer(intg), parameter part_deriv_s2_s3_s4
Cross derivative in the s2, s3 and s4 direction i.e., d^3u/ds2ds3ds4.
subroutine basis_quadrature_order_set_ptr(BASIS, ORDER, ERR, ERROR,)
Sets/changes the order of a quadrature for a basis quadrature identified by a pointer.
integer(intg), parameter part_deriv_s1_s3
Cross derivative in the s1 and s3 direction i.e., d^2u/ds1ds3.
subroutine, public basis_quadrature_number_of_gauss_xi_get(BASIS, QUADRATURE_NUMBER_OF_GAUSS_XI, ERR, ERROR,)
Get the number of Gauss points in each xi direction on a basis quadrature identified by a pointer...
integer(intg), parameter, public basis_multiquartic_radial_interpolation
Multiquartic Radial interpolation specification.
integer(intg), parameter third_part_deriv
Third partial derivative i.e., d^3u/ds^3.
subroutine basis_lhtp_family_create(BASIS, ERR, ERROR,)
Creates and initialises a Lagrange-Hermite tensor product basis family that has already been allocate...
A buffer type to allow for an array of pointers to a BASIS_TYPE.
subroutine basis_number_of_xi_set_number(USER_NUMBER, NUMBER_OF_XI, ERR, ERROR,)
Sets/changes the number of xi directions where the basis is identified by user number.
subroutine basis_quadrature_create(BASIS, ERR, ERROR,)
Creates the quadrature and quadrature schemes on a basis.
real(dp) function hermite_quadratic_evaluate(NODE_INDEX, NODE_DERIVATIVE_INDEX, PARTIAL_DERIVATIVE_INDEX, SPECIAL_NODE_INDEX, XI, ERR, ERROR)
Evaluates a 1D quadratic Hermite basis function.
integer(intg), parameter, public basis_lagrange_hermite_tp_type
Lagrange-Hermite tensor product basis type.
integer(intg), parameter, public basis_gauss_simplex_quadrature
Gauss-Legendre for Simplex elements quadrature.
logical, save, public diagnostics1
.TRUE. if level 1 diagnostic output is active in the current routine
integer(intg), parameter part_deriv_s1_s4
Cross derivative in the s1 and s4 direction i.e., d^2u/ds1ds4.
A buffer type to allow for an array of pointers to a QUADRATURE_SCHEME_TYPE.
subroutine, public basis_type_get(BASIS, TYPE, ERR, ERROR,)
get the type for a basis is identified by a a pointer.
integer(intg), parameter, public basis_collapsed_at_xi1
The Xi direction at the xi=1 end of this Xi direction is collapsed.
subroutine basis_sub_basis_create(PARENT_BASIS, NUMBER_OF_XI, XI_DIRECTIONS, SUB_BASIS, ERR, ERROR,)
Creates a sub-basis on a parent basis.
Contains information for a particular quadrature scheme.
integer(intg), parameter, public basis_cubic_lagrange_interpolation
Cubic Lagrange interpolation specification.
subroutine basis_collapsed_xi_set_ptr(BASIS, COLLAPSED_XI, ERR, ERROR,)
Sets/changes the collapsed xi flags for a basis is identified by a a pointer.
subroutine basis_simplex_basis_create(BASIS, ERR, ERROR,)
Creates and initialises a simplex basis that has already been allocated BASIS_CREATE_START.
subroutine, public basis_quadrature_type_get(BASIS, QUADRATURE_TYPE, ERR, ERROR,)
get the quadrature type on a basis identified by a pointer.
real(dp), parameter twopi
The double value of 2pi.
Interpolates the requested partial derivative index(ices) of the element parameters for basis functio...
integer(intg), parameter part_deriv_s2_s4
Cross derivative in the s2 and s4 direction i.e., d^2u/ds2ds4.
Evaluates the list of gauss points and weights for a given basis type and order.
integer(intg), parameter, public basis_auxilliary_type
Auxillary basis type.
integer(intg), parameter, public diagnostic_output_type
Diagnostic output type.
subroutine gauss_legendre(N, ALPHA, BETA, X, W, ERR, ERROR,)
This routine calculates the weights and abscissae for a Gauss-Legendre quadrature scheme...
integer(intg), parameter, public basis_quadratic1_hermite_interpolation
Quadratic Hermite (no derivative at xi=0) interpolation specification.
recursive subroutine basis_family_number_find(USER_NUMBER, FAMILY_NUMBER, BASIS, ERR, ERROR,)
Finds and returns in BASIS a pointer to the basis with the given USER_NUMBER and FAMILY_NUMBER. If no basis with that number and family number exists then BASIS is returned nullified.
integer(intg), parameter part_deriv_s1_s2
Cross derivative in the s1 and s2 direction i.e., d^2u/ds1ds2.
integer(intg), parameter part_deriv_s2_s2
Second partial derivative in the s2 direction i.e., d^2u/ds2ds2.
integer(intg), parameter, public basis_low_quadrature_scheme
Identifier for a low order quadrature scheme.
subroutine, public basis_number_of_xi_get(BASIS, NUMBER_OF_XI, ERR, ERROR,)
Gets the number of xi directions for a basis.
subroutine basis_interpolation_xi_set_ptr(BASIS, INTERPOLATION_XI, ERR, ERROR,)
Sets/changes the interpolation type in each xi directions for a basis identified by a pointer...
integer(intg), parameter part_deriv_s3_s3
Second partial derivative in the s3 direction i.e., d^2u/ds3ds3.
integer(intg), parameter basis_fourier_interpolation
Fourier interpolation.
real(dp) function basis_interpolate_xi_dp(BASIS, PARTIAL_DERIV_INDEX, XI, ELEMENT_PARAMETERS, ERR, ERROR)
Interpolates the appropriate partial derivative index of the element parameters at position XI for th...
real(dp) function basis_lhtp_basis_evaluate_dp(BASIS, NODE_NUMBER, DERIVATIVE_NUMBER, PARTIAL_DERIV_INDEX, XI, ERR, ERROR)
Evaluates the double precision Lagrange/Hermite/Fourier tensor product basis function for the given B...
real(dp) function lagrange_quadratic_evaluate(NODE_INDEX, PARTIAL_DERIVATIVE_INDEX, XI, ERR, ERROR)
Evaluates a 1D quadratic Lagrange basis function.
real(dp) function lagrange_cubic_evaluate(NODE_INDEX, PARTIAL_DERIVATIVE_INDEX, XI, ERR, ERROR)
Evaluates a 1D cubic Lagrange basis function.
subroutine basis_quadrature_order_set_number(USER_NUMBER, ORDER, ERR, ERROR,)
Sets/changes the order of a quadrature for a basis quadrature identified by a user number...
integer(intg), parameter, public basis_cubic_hermite_interpolation
Cubic Hermite interpolation specification.
real(dp) function hermite_cubic_evaluate(NODE_INDEX, NODE_DERIVATIVE_INDEX, PARTIAL_DERIVATIVE_INDEX, XI, ERR, ERROR)
Evaluates a 1D cubic Hermite basis function.
integer(intg), parameter part_deriv_s2_s4_s4
Cross derivative in the s2, s4 and s4 direction i.e., d^3u/ds2ds4^2.
Contains all information about a basis .
Sets/changes the order of a quadrature for a basis quadrature.
Contains information on the defined basis functions.
integer(intg), parameter, public basis_quadratic2_interpolation_order
Quadratic (no derivative at xi=1) interpolation order.
Sets/changes the type for a basis.
Sets/changes the quadrature type for a basis.
subroutine, public basis_create_start(USER_NUMBER, BASIS, ERR, ERROR,)
Starts the creation of a new basis The default values of the BASIS attributes are: ...
subroutine basis_quadrature_initialise(BASIS, ERR, ERROR,)
Initialises a quadrature on the given basis.
Flags an error condition.
integer(intg), dimension(23, 4) partial_derivative_index
Partial derivative index map. PARTIAL_DERIVATIVE_INDEX(idx,nic) gives the order of the partial deriva...
integer(intg), parameter, public basis_linear_lagrange_interpolation
Linear Lagrange interpolation specification.
subroutine, public bases_finalise(ERR, ERROR,)
Finalises the bases and deallocates all memory.
real(dp) function lagrange_linear_evaluate(NODE_INDEX, PARTIAL_DERIVATIVE_INDEX, XI, ERR, ERROR)
Evaluates a 1D linear Lagrange basis function.
subroutine basis_quadrature_type_set_ptr(BASIS, TYPE, ERR, ERROR,)
Sets/changes the quadrature type on a basis identified by a pointer.
subroutine, public basis_quadrature_order_get(BASIS, QUADRATURE_ORDER, ERR, ERROR,)
Get the order of a quadrature for a basis quadrature identified by a pointer.
integer(intg), parameter, public basis_linear_interpolation_order
Linear interpolation order.
This module contains all kind definitions.
subroutine, public basis_local_node_xi_calculate(BASIS, LOCAL_NODE_NUMBER, XI, ERR, ERROR,)
Calculates the xi location of a local node in a basis.
integer(intg), parameter part_deriv_s1_s2_s3
Cross derivative in the s1, s2 and s3 direction i.e., d^3u/ds1ds2ds3.
integer(intg), parameter part_deriv_s1_s2_s4
Cross derivative in the s1, s2 and s4 direction i.e., d^3u/ds1ds2ds4.
integer(intg), parameter basis_hermite_interpolation
Hermite interpolation.