From 3b1555ffc32fe374545c1e886413871b37486768 Mon Sep 17 00:00:00 2001 From: "Hrant P. Hratchian" Date: Thu, 20 Jun 2024 18:20:21 -0700 Subject: [PATCH] Updates to a few modules. --- src/mqc_algebra2.F03 | 90 ++++++++++++++++++++++++++++++++++++++++-- src/mqc_gaussian.F03 | 49 +++++++++++++++++++++++ src/mqc_general.F03 | 25 ++++++++++++ src/mqc_matwrapper.F03 | 2 +- 4 files changed, 161 insertions(+), 5 deletions(-) diff --git a/src/mqc_algebra2.F03 b/src/mqc_algebra2.F03 index e40819bb..a94009c5 100644 --- a/src/mqc_algebra2.F03 +++ b/src/mqc_algebra2.F03 @@ -129,6 +129,11 @@ Module MQC_Algebra2 module procedure MQC_Variable_dot_product end interface ! +! Interface MQC Varaiable vector-vector outer product --> outer_product. + interface outer_product + module procedure MQC_Variable_outer_product + end interface +! ! Interface MQC Varaiable full SUM. interface sum module procedure MQC_Variable_Array_Sum @@ -3344,10 +3349,6 @@ function MQC_Variable_dot_product(mqcVariable1,mqcVariable2) result(mqcVariableO call mqc_error_i('dot_product must be between two vectors.',6, & 'RANK(mqcVariable1)',mqcVariable1%getRank(), & 'RANK(mqcVariable2)',mqcVariable2%getRank()) - if(.not.MQC_Variable_isConformable(mqcVariable1,mqcVariable2)) & - call mqc_error_l('MQC Division only allowed for conformable values.',6, & - 'MQC_Variable_isConformable(mqcVariable1,mqcVariable2)', & - MQC_Variable_isConformable(mqcVariable1,mqcVariable2)) ! ! Calculate the dot product. ! @@ -3357,6 +3358,87 @@ function MQC_Variable_dot_product(mqcVariable1,mqcVariable2) result(mqcVariableO end function MQC_Variable_dot_product +! +!PROCEDURE MQC_Variable_outer_product + function MQC_Variable_outer_product(mqcVariable1,mqcVariable2) result(mqcVariableOut) +! +! This function solves the outer product of two vectors stored in the input +! dummy MQC_Variables arguments and . If these +! two arguments do not correspond to vectors, this function fails. +! +! +! H. P. Hratchian, 2024. +! +! +! Variable Declarations. + class(MQC_Variable),intent(in)::mqcVariable1,mqcVariable2 + type(MQC_Variable)::mqcVariableOut + integer(kind=int64)::typeCode1,typeCode2,len1,len2 + integer(kind=int64),dimension(:),allocatable::vectorInteger1, & + vectorInteger2 + integer(kind=int64),dimension(:,:),allocatable::matrixInteger + real(kind=real64),dimension(:),allocatable::vectorReal1, & + vectorReal2 + real(kind=real64),dimension(:,:),allocatable::matrixReal +! +! +! Start by ensuring the two MQC Variables are vectors and are conformable. +! + if(mqcVariable1%getRank().ne.1.or.mqcVariable2%getRank().ne.1) & + call mqc_error_i('outer_product must be between two vectors.',6, & + 'RANK(mqcVariable1)',mqcVariable1%getRank(), & + 'RANK(mqcVariable2)',mqcVariable2%getRank()) +! +! Determine the type codes for the two vectors and ensure they're allowable. +! + typeCode1 = MQC_Variable_getTypeCode(mqcVariable1) + typeCode2 = MQC_Variable_getTypeCode(mqcVariable2) + if(typeCode1.eq.0) & + call mqc_error_I('Outer product: Var1 is of UNKONWN type.', 6, & + 'typeCode1', typeCode1) + if(typeCode2.eq.0) & + call mqc_error_I('Outer product: Var2 is of UNKONWN type.', 6, & + 'typeCode2', typeCode2) +! +! Now, do the outer product and set the output function value accordingly. +! + len1 = SIZE(mqcVariable1) + len2 = SIZE(mqcVariable2) + select case(typeCode1*10 + typeCode2) + case(22) + vectorReal1 = mqcVariable1 + vectorReal2 = mqcVariable2 + matrixReal = matmul(RESHAPE(vectorReal1,[ len1,1 ]), & + RESHAPE(vectorReal2,[ 1,len2 ])) + mqcVariableOut = matrixReal + case(23) + vectorReal1 = mqcVariable1 + vectorInteger2 = mqcVariable2 + matrixReal = matmul(RESHAPE(vectorReal1,[ len1,1 ]), & + RESHAPE(vectorInteger2,[ 1,len2 ])) + mqcVariableOut = matrixReal + case(32) + vectorInteger1 = mqcVariable1 + vectorReal2 = mqcVariable2 + matrixReal = matmul(RESHAPE(vectorInteger1,[ len1,1 ]), & + RESHAPE(vectorReal2,[ 1,len2 ])) + mqcVariableOut = matrixReal + case(33) + vectorInteger1 = mqcVariable1 + vectorInteger2 = mqcVariable2 + matrixInteger = matmul(RESHAPE(vectorInteger1,[ len1,1 ]), & + RESHAPE(vectorInteger2,[ 1,len2 ])) + mqcVariableOut = matrixInteger + case default + call mqc_error_I('Outer product: Combination of var1 and var 2 types is UNKNOWN.', 6, & + 'typeCode1', typeCode1, & + 'typeCode2', typeCode2 ) + end select +! + return + end function MQC_Variable_outer_product + + ! !PROCEDURE MQC_Variable_MatrixVector function MQC_Variable_MatrixVector(mqcMatrix,mqcVector) result(mqcVariableOut) diff --git a/src/mqc_gaussian.F03 b/src/mqc_gaussian.F03 index 178a3b38..f1eafdbb 100644 --- a/src/mqc_gaussian.F03 +++ b/src/mqc_gaussian.F03 @@ -116,6 +116,11 @@ Module MQC_Gaussian !---------------------------------------------------------------- ! ! +! Interface BasisType2FunctionInfo... + interface basisType2FunctionInfo + module procedure MQC_Gaussian_Unformatted_Matrix_Basis_Type_To_Function_Info + end interface +! ! !---------------------------------------------------------------- ! | @@ -2458,6 +2463,50 @@ Function MQC_Gaussian_Unformatted_Matrix_Get_Basis_Info_Array(fileinfo,label) & return end Function MQC_Gaussian_Unformatted_Matrix_Get_Basis_Info_Array + +!===================================================================== +! +!PROCEDURE MQC_Gaussian_Unformatted_Matrix_Basis_Type_To_Function_Info + Subroutine MQC_Gaussian_Unformatted_Matrix_Basis_Type_To_Function_Info(basisType, & + angularMomentum,component,isCartesian) +! +! This routine is used to interpret an element of the basis function type +! list. The basis type value is the only input argument (see functions +! MQC_Gaussian_Unformatted_Matrix_Get_Basis_Info_Element and +! MQC_Gaussian_Unformatted_Matrix_Get_Basis_Info_Array to get the basis type +! value for a basis function) The output of this routine is the angular +! momentum of the basis function, the component (x, y, z for l=1, for +! example), and a logical that is TRUE for Cartesian basis functions and +! FALSE for pure basis functions. +! +! +! H. P. Hratchian, 2024. +! +! +! Variable Declarations. +! + implicit none + + integer(kind=int64),intent(in)::basisType + integer,intent(out),optional::angularMomentum,component + logical,intent(out),optional::isCartesian + integer::myBasisType,myAngularMomentum,myComponent +! +! Decode the basisType number... +! + myBasisType = ABS(basisType) + myAngularMomentum = myBasisType/1000 + myComponent = mod(myBasisType,1000) +! + if(PRESENT(angularMomentum)) angularMomentum = myAngularMomentum + if(PRESENT(component)) component = myComponent + if(PRESENT(isCartesian)) isCartesian = basisType.gt.0 +! + return + end subroutine MQC_Gaussian_Unformatted_Matrix_Basis_Type_To_Function_Info + + + !===================================================================== ! diff --git a/src/mqc_general.F03 b/src/mqc_general.F03 index aa046f80..2c2d3111 100644 --- a/src/mqc_general.F03 +++ b/src/mqc_general.F03 @@ -663,6 +663,31 @@ Subroutine mqc_get_command_argument_integer(argNum,iArgument) Return End Subroutine mqc_get_command_argument_integer +! +!PROCEDURE mqc_get_command_argument_real + Subroutine mqc_get_command_argument_real(argNum,rArgument) +! +! This subroutine is used to get a real argument from the command line. The +! output dummy argument rArgument is returned with the real number found as +! command line arugment argNum. If this command line argument is NOT an +! real, the program will SegFault. +! +! -H. P. Hratchian, 2024. +! +! + implicit none + integer(kind=int64),intent(in)::argNum + real(kind=real64),intent(out)::rArgument + character(len=:),allocatable::argument +! +! Do the work... +! + call mqc_get_command_argument(argNum,argument) + read(argument,*) rArgument +! + Return + End Subroutine mqc_get_command_argument_real + ! !---------------------------------------------------------------- diff --git a/src/mqc_matwrapper.F03 b/src/mqc_matwrapper.F03 index 6ce47d57..8e1443f8 100644 --- a/src/mqc_matwrapper.F03 +++ b/src/mqc_matwrapper.F03 @@ -1,4 +1,4 @@ -!hph #define UCMGAUOPEN +!hph#define UCMGAUOPEN Module MQC_MatWrapper ! ! **********************************************************************