Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/prop-base/pde_solver_control.f90.svn-base

+K 13
+svn:eol-style
+V 6
+native
+END

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/prop-base/pde_solver_control_test.f90.svn-base

+K 13
+svn:eol-style
+V 6
+native
+END

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/prop-base/system_constants.f90.svn-base

+K 13
+svn:eol-style
+V 6
+native
+END

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/prop-base/test_driver.f90.svn-base

+K 13
+svn:eol-style
+V 6
+native
+END

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/prop-base/test_straight.xlsx.svn-base

+K 13
+svn:mime-type
+V 24
+application/octet-stream
+END

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/prop-base/test_straight_coarse.elem.svn-base

+K 13
+svn:eol-style
+V 6
+native
+END

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/prop-base/test_straight_coarse.node.svn-base

+K 13
+svn:eol-style
+V 6
+native
+END

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/prop-base/test_straight_coarse.time.svn-base

+K 13
+svn:eol-style
+V 6
+native
+END

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/prop-base/vector_def.f90.svn-base

+K 13
+svn:eol-style
+V 6
+native
+END

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/prop-base/vector_test.f90.svn-base

+K 13
+svn:eol-style
+V 6
+native
+END

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/text-base/Makefile.svn-base

+# Fortran compiler
+FF = gfortran
+
+# general compiler flags
+FLAGS =
+
+# compiler flags for linear algebra
+LINALG = -llapack -lblas
+
+# list of required object files
+OBJECTS = fruit_util.o fruit.o \
+		system_constants.o \
+		log_messages.o log_message_control.o \
+		vector_def.o dense_matrix_def.o band_sym_matrix_def.o \
+        linear_solver.o \
+        material_data.o field_data.o boundary_data.o \
+		pde_solver_constants.o \
+		constitutive.o kinematic.o \
+		material_model.o \
+		interpolation.o integration.o \
+		pde_solver_control.o \
+		log_messages_test.o log_message_control_test.o \
+		vector_test.o dense_matrix_test.o band_sym_matrix_test.o \
+        linear_solver_test.o \
+		material_data_test.o field_data_test.o boundary_data_test.o \
+		constitutive_test.o \
+		pde_solver_control_test.o \
+		test_driver.o
+
+all: test_driver check
+
+check:
+	./test_driver.a
+	
+doc:
+	doxygen doc-dynsws
+
+test_driver: $(OBJECTS)
+	$(FF) -o test_driver.a $(OBJECTS) $(LINALG) $(FLAGS)
+
+		
+		
+fruit_util.o: fruit_util.f90
+	$(FF) -c fruit_util.f90 $(FLAGS)
+	
+fruit.o: fruit.f90
+	$(FF) -c fruit.f90 $(FLAGS)
+	
+	
+	
+test_driver.o: test_driver.f90
+	$(FF) -c test_driver.f90 $(LINALG) $(FLAGS)
+	
+	
+	
+band_sym_matrix_def.o: band_sym_matrix_def.f90
+	$(FF) -c band_sym_matrix_def.f90 $(FLAGS)
+	
+band_sym_matrix_test.o: band_sym_matrix_test.f90
+	$(FF) -c band_sym_matrix_test.f90 $(FLAGS)
+	
+	
+	
+boundary_data.o: boundary_data.f90
+	$(FF) -c boundary_data.f90 $(FLAGS)
+	
+boundary_data_test.o: boundary_data_test.f90
+	$(FF) -c boundary_data_test.f90 $(FLAGS)
+	
+
+	
+constitutive.o: constitutive.f90
+	$(FF) -c constitutive.f90 $(FLAGS)
+	
+constitutive_test.o: constitutive_test.f90
+	$(FF) -c constitutive_test.f90 $(FLAGS)
+	
+	
+	
+dense_matrix_def.o: dense_matrix_def.f90
+	$(FF) -c dense_matrix_def.f90 $(FLAGS)
+	
+dense_matrix_test.o: dense_matrix_test.f90
+	$(FF) -c dense_matrix_test.f90 $(FLAGS)
+	
+	
+	
+field_data.o: field_data.f90
+	$(FF) -c field_data.f90 $(FLAGS)
+	
+field_data_test.o: field_data_test.f90
+	$(FF) -c field_data_test.f90 $(FLAGS)
+	
+	
+	
+integration.o: integration.f90
+	$(FF) -c integration.f90 $(FLAGS)
+	
+	
+	
+interpolation.o: interpolation.f90
+	$(FF) -c interpolation.f90 $(FLAGS)
+	
+	
+	
+kinematic.o: kinematic.f90
+	$(FF) -c kinematic.f90 $(FLAGS)
+	
+	
+    
+linear_solver.o: linear_solver.f90
+	$(FF) -c linear_solver.f90 $(LINALG) $(FLAGS)
+	
+linear_solver_test.o: linear_solver_test.f90
+	$(FF) -c linear_solver_test.f90 $(LINALG) $(FLAGS)
+	
+    
+	
+log_message_control.o: log_message_control.f90
+	$(FF) -c log_message_control.f90 $(FLAGS)
+	
+log_message_control_test.o: log_message_control_test.f90
+	$(FF) -c log_message_control_test.f90 $(FLAGS)
+	
+	
+	
+log_messages.o: log_messages.f90
+	$(FF) -c log_messages.f90 $(FLAGS)
+	
+log_messages_test.o: log_messages_test.f90
+	$(FF) -c log_messages_test.f90 $(FLAGS)
+	
+	
+	
+material_model.o: material_model.f90
+	$(FF) -c material_model.f90 $(FLAGS)
+	
+	
+	
+material_data.o: material_data.f90
+	$(FF) -c material_data.f90 $(FLAGS)
+	
+material_data_test.o: material_data_test.f90
+	$(FF) -c material_data_test.f90 $(FLAGS)
+	
+	
+	
+pde_solver_constants.o: pde_solver_constants.f90
+	$(FF) -c pde_solver_constants.f90 $(FLAGS)
+	
+	
+	
+pde_solver_control.o: pde_solver_control.f90
+	$(FF) -c pde_solver_control.f90 $(FLAGS)
+	
+pde_solver_control_test.o: pde_solver_control_test.f90
+	$(FF) -c pde_solver_control_test.f90 $(FLAGS)
+	
+	
+	
+system_constants.o: system_constants.f90
+	$(FF) -c system_constants.f90 $(FLAGS)
+	
+	
+	
+vector_def.o: vector_def.f90
+	$(FF) -c vector_def.f90 $(FLAGS)
+	
+vector_test.o: vector_test.f90
+	$(FF) -c vector_test.f90 $(FLAGS)
+	
+
+
+clean:
+	rm -rf *.o
+	
+clean_test:
+	rm -rf test_driver.o
+	
+clean_boundary_data: 	clean_test \
+						clean_integration clean_interpolation \
+						clean_pde_solver_control
+	rm -rf boundary_data.o boundary_data_test.o
+	
+clean_bsm: 	clean_test \
+			clean_pde_solver_control
+	rm -rf band_sym_matrix_def.o band_sym_matrix_test.o
+	
+clean_constitutive: clean_test \
+					clean_integration clean_material_model
+	rm -rf constitutive.o constitutive_test.o
+	
+clean_dm: 	clean_test \
+			clean_bsm clean_constitutive clean_integration \
+			clean_interpolation clean_kinematic clean_material_model \
+			clean_pde_solver_control
+	rm -rf dense_matrix_def.o dense_matrix_test.o
+	
+clean_field_data: 	clean_test \
+					clean_integration clean_interpolation \
+					clean_kinematic clean_boundary_data \
+					clean_pde_solver_control
+	rm -rf field_data.o field_data_test.o
+	
+clean_integration: 	clean_test \
+					clean_pde_solver_control
+	rm -rf integration.o
+	
+clean_interpolation: 	clean_test \
+						clean_integration
+	rm -rf interpolation.o
+	
+clean_kinematic: 	clean_test \
+					clean_integration
+	rm -rf kinematic.o
+	
+clean_linear_solver: clean_test
+	rm -rf linear_solver.o linear_solver_test.o
+	
+clean_log: 	clean_test \
+			clean_bsm clean_dm clean_vec \
+			clean_material_data clean_field_data \
+			clean_constitutive \
+			clean_pde_solver_control
+	rm -rf log_message_control.o log_message_control_test.o
+	
+clean_material_data: 	clean_test \
+						clean_integration clean_field_data \
+						clean_pde_solver_control
+	rm -rf material_data.o material_data_test.o
+	
+clean_material_model: clean_test
+	rm -rf material_model.o
+	
+clean_msg: 	clean_test \
+			clean_bsm clean_dm clean_log clean_vec \
+			clean_material_data clean_field_data \
+			clean_constitutive \
+			clean_pde_solver_control
+	rm -rf log_messages.o log_messages_test.o
+	
+clean_pde_solver_constants: clean_test \
+							clean_integration \
+							clean_pde_solver_control
+	rm -rf pde_solver_constants.o
+	
+clean_pde_solver_control: 	clean_test
+	rm -rf pde_solver_control.o pde_solver_control_test.o
+	
+clean_system_constants: clean_test \
+						clean_log clean_msg \
+						clean_constitutive clean_kinematic \
+						clean_interpolation \
+						clean_pde_solver_constants \
+						clean_material_data clean_field_data \
+						clean_pde_solver_control
+	rm -rf system_constants.o
+	
+clean_vec: 	clean_test \
+			clean_bsm clean_dm clean_integration \
+			clean_material_model \
+			clean_pde_solver_control
+	rm -rf vector_def.o vector_test.o

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/text-base/README.txt.svn-base

+* ------------------------------------------------------------------------
+*
+* DynSWS 1.0
+* Dynamic Model of Soil-Water-Structure Interaction
+*
+* Contact:
+* Brandon Karchewski
+* c/o Dept. of Civil Engineering (JHE-301)
+* McMaster University
+* 1280 Main Street West
+* Hamilton, Ontario, Canada
+* L8S 4L7
+* T: (905) 529-6569
+* E: karcheba@mcmaster.ca
+*
+* ------------------------------------------------------------------------
+
+
+------------------
+REVISION HISTORY:
+------------------
+
+    [26 March 2012]
+        -Initialized this README file
+    
+    [29 March 2012]
+        -Added build instructions for prerequisite libraries
+        -Modified build instructions for DynSWS
+    
+    [1 April 2012]
+        -Added descriptions of new modules
+        
+    [2 April 2012]
+        -Added descriptions of new modules
+    
+    [3 April 2012]
+        -Added descriptions of new modules
+        -Modified build instructions for clarity
+        
+    [5 April 2012]
+        -Added documentation instructions
+        
+    [6 April 2012]
+        -Added descriptions of new modules
+        
+        
+        
+------------------------
+BACKGROUND INFORMATION:
+------------------------
+
+This file relates to the program DynSWS. New users/developers should see
+the following documents to become acquainted with the project:
+
+Software Requirements Specification
+DynSWS-SRS-1.0.pdf
+
+Module Guide:
+DynSWS-MG-1.0.pdf
+
+Module Interface Specification:
+DynSWS-MIS-1.0.pdf
+
+SRS, MG, and MIS for PDE Solver Module:
+DynSWS-PDE-1.0.pdf
+
+
+
+-----------------------
+PROJECT PREREQUISITES:
+-----------------------
+
+The following netlib linear algebra packages should be available on
+the build path to be linked against:
+
+    BLAS (Basic Linear Algebra Subprograms)
+        http://www.netlib.org/blas/
+
+    LAPACK (Linear Algebra PACKage)
+        http://www.netlib.org/lapack/
+
+Both can be obtained by visiting the LAPACK project website and
+downloading the latest version of LAPACK. As of this moment it
+is LAPACK 3.4.0. Here is a minimal instruction set for building
+the package (assuming the compiler is gfortran):
+
+    1. Extract downloaded package to desired directory.
+    2. cd to the root directory of the LAPACK package.
+    3. Delete the .example extension from make.inc.example.
+        -Note: The user may want to look inside this file to make
+                sure the configuration is appropriate for their
+                system.
+    4. cd BLAS/SRC
+    5. make
+    6. cd ../../
+    7. make
+
+This should build and test the reference BLAS and the LAPACK
+package. The user may want to use their own machine-tuned version,
+but that will be left at their discretion. At any rate, the
+instructions above should leave the user with a generic BLAS and
+LAPACK build. The libraries generated will have the filenames:
+    
+    librefblas.a
+    libtmglib.a
+    liblapack.a
+    
+The user should copy these into a folder that is on the build
+path for their compiler (generally a folder called 'lib' within
+the file system of the shell program or the compiler suite).
+Generally, the user should also change the name of 'librefblas.a'
+to 'libblas.a' so that the compiler can find it with the -lblas
+switch.
+       
+       
+       
+-----------------------------
+DOCUMENTATION PREREQUISITES:
+-----------------------------
+
+To build the documentation from the source code, the doxygen system is
+required. This can be obtained from:
+    
+    http://www.stack.nl/~dimitri/doxygen/
+    
+The user should follow the installation instructions associated with their
+system. Note that this item is not required to build the actual program.
+
+
+
+------------------
+PROJECT CONTENTS:
+------------------
+
+    README.txt
+        Contains general information about project and build instructions
+
+    Makefile
+        Makefile for building the program including unit tests
+
+    fruit.f90
+    fruit_util.f90
+        Unit testing framework (FRUIT) for Fortran
+        See http://sourceforge.net/apps/mediawiki/fortranxunit/index.php?title=Main_Page
+            http://sourceforge.net/projects/fortranxunit/
+            http://fortranxunit.sourceforge.net/
+
+    band_sym_matrix_def.f90
+    band_sym_matrix_test.f90
+        Banded Symmetric Matrix ADT and associated unit tests
+        
+    boundary_data.f90
+    boundary_data_test.f90
+        Boundary Data Module and associated tests
+        
+    constitutive.f90
+        Constitutive Matrix module
+
+    dense_matrix_def.f90
+    dense_matrix_test.f90
+        Dense Matrix ADT and associated unit tests
+        
+    field_data.f90
+    field_data_test.f90
+        Field Data Module and associated tests
+        
+    integration.f90
+        Body Element Integration and Traction Element Integration modules
+        
+    interpolation.f90
+        Body Element Interpolation and Traction Element Interpolation modules
+        
+    kinematic.f90
+        Kinematic Matrix Module
+
+    log_message_control.f90
+    log_message_control_test.f90
+        Log Message Control module and associated unit tests
+
+    log_messages.f90
+    log_messages_test.f90
+        Log Messages Module and associated unit tests
+        
+    material_data.f90
+    material_data_test.f90
+        Material Property Data module and associated unit tests
+        
+    pde_solver_control.f90
+    pde_solver_control_test.f90
+        PDE Solver Control module and associated unit tests
+
+    system_constants.f90
+        System Constants Module
+
+    test_driver.f90
+        Driver program for unit tests
+
+    vector_def.f90
+    vector_test.f90
+        Vector ADT and associated unit tests
+        
+        
+        
+--------------------
+BUILD INSTRUCTIONS:
+--------------------
+
+To build and run the test suite:
+    1. Open shell
+    2. cd to directory of DynSWS
+    3. Open Makefile and set FF to the name of your Fortran compiler
+            (the default is gfortran)
+    4. make
+    
+This will build the program and run the unit test suite. The results of
+the unit tests will be printed to the console. The test suite will also
+generate testName.log and test_straight_coarse.log files. Assuming that
+the unit tests all pass, the contents are not important as these are just
+by-products of the exception checking tests.
+
+
+
+--------------------
+BUILD DOCUMENTATION:
+--------------------
+
+To build documentation from the source code:
+    1. Open shell
+    2. cd to directory of DynSWS
+    3. make doc
+    
+Note that this generates both html and latex format documentation. The
+html version will work out of the box and is contained in the directory
+<dir of DynSWS>/doc/html. The latex documentation may be build as follows
+(assuming the shell is still open and the current directory is that of
+DynSWS):
+
+    1. cd doc/latex
+    2. make
+    
+Note that the latex version has issues at this point in time. The author's
+recommended version of the documentation is the html version.

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/text-base/band_sym_matrix_def.f90.svn-base

+! ------------------------------------------------------------------------
+!>  \brief Module defining Banded Symmetric Matrix data type
+! ------------------------------------------------------------------------
+MODULE band_sym_matrix_def
+USE log_message_control   !> Printing log/error messages
+USE log_messages          !> Log/error codes
+USE vector_def            !> Vector data type
+USE dense_matrix_def      !> Dense Matrix data type
+
+IMPLICIT NONE
+PRIVATE
+
+
+
+! ************************************************************************
+! EXPORTS
+! ************************************************************************
+
+
+!> Exported data types
+PUBLIC :: bandSymMatrixT
+
+!> Exported interfaces
+PUBLIC :: bsm_init, bsm_clean, &
+          bsm_numRows, bsm_halfBW, &
+          bsm_get, bsm_set, &
+          bsm_setDecomp, bsm_isDecomposed, &
+          bsm_add, OPERATOR (+), bsm_mappedAdd, &
+          bsm_scalMul, OPERATOR (*), bsm_vecMul
+
+
+
+! ************************************************************************
+! LOCAL CONSTANTS
+! ************************************************************************
+
+!> Sender code for Banded SYmmetric MATrix module
+INTEGER, PARAMETER :: sdr = BSYMAT
+
+
+
+! ************************************************************************
+! DATA TYPES
+! ************************************************************************
+
+
+! ------------------------------------------------------------------------
+!>  \brief Banded Symmetric Matrix ADT structure
+!!
+!!  \param dat            Array to store packed (banded symmetric) data
+!!  \param decomp         Array to store packed Cholesky decomposition
+!!  \param is_decomposed  Switch indicating whether decomposition is present
+! ------------------------------------------------------------------------
+TYPE bandSymMatrixT
+  DOUBLE PRECISION, ALLOCATABLE :: dat(:,:)
+  DOUBLE PRECISION, ALLOCATABLE :: decomp(:,:)
+  LOGICAL :: is_decomposed
+END TYPE bandSymMatrixT
+
+
+
+! ************************************************************************
+! INTERFACES
+! ************************************************************************
+
+
+!>  \brief Interface to constructor
+INTERFACE bsm_init
+  MODULE PROCEDURE bsm_init_
+  MODULE PROCEDURE bsm_init_exc_
+END INTERFACE bsm_init
+
+!>  \brief Interface to destructor
+INTERFACE bsm_clean
+  MODULE PROCEDURE bsm_clean_
+END INTERFACE bsm_clean
+
+!>  \brief Interface to number of rows
+INTERFACE bsm_numRows
+  MODULE PROCEDURE bsm_num_rows_
+END INTERFACE bsm_numRows
+
+!>  \brief Interface to half bandwidth
+INTERFACE bsm_halfBW
+  MODULE PROCEDURE bsm_half_bw_
+END INTERFACE bsm_halfBW
+
+!>  \brief Interface to getter for individual entries in matrix
+INTERFACE bsm_get
+  MODULE PROCEDURE bsm_get_
+  MODULE PROCEDURE bsm_get_exc_
+END INTERFACE bsm_get
+
+!>  \brief Interface to setter for individual entries in matrix
+INTERFACE bsm_set
+  MODULE PROCEDURE bsm_set_
+  MODULE PROCEDURE bsm_set_exc_
+END INTERFACE bsm_set
+
+!>  \brief Interface to setter for Cholesky decomposition
+INTERFACE bsm_setDecomp
+  MODULE PROCEDURE bsm_set_decomp_
+  MODULE PROCEDURE bsm_set_decomp_exc_
+END INTERFACE bsm_setDecomp
+
+!>  \brief Interface to check if decomposition exists
+INTERFACE bsm_isDecomposed
+  MODULE PROCEDURE bsm_is_decomposed_
+END INTERFACE bsm_isDecomposed
+
+!>  \brief Interface for addition of two bandSymMatrixT
+INTERFACE bsm_add
+  MODULE PROCEDURE bsm_add_
+  MODULE PROCEDURE bsm_add_exc_
+END INTERFACE bsm_add
+
+!>  \brief Operator overload for addition
+INTERFACE OPERATOR (+)
+  MODULE PROCEDURE bsm_add_
+END INTERFACE
+
+!>  \brief Interface for adding a matrixT to a bandSymMatrixT in place with mapping indices
+INTERFACE bsm_mappedAdd
+  MODULE PROCEDURE bsm_mapped_add_
+  MODULE PROCEDURE bsm_mapped_add_exc_
+END INTERFACE bsm_mappedAdd
+
+!>  \brief Interface for scalar multiplication
+INTERFACE bsm_scalMul
+  MODULE PROCEDURE bsm_scal_mul_mat_scal_
+END INTERFACE bsm_scalMul
+
+!>  \brief Interface for matrix-vector multiplication
+INTERFACE bsm_vecMul
+  MODULE PROCEDURE bsm_vec_mul_
+  MODULE PROCEDURE bsm_vec_mul_exc_
+END INTERFACE bsm_vecMul
+
+!>  \brief Operator overload for multiplication (includes scalar and vector multiplication)
+INTERFACE OPERATOR (*)
+  MODULE PROCEDURE bsm_scal_mul_mat_scal_
+  MODULE PROCEDURE bsm_scal_mul_scal_mat_
+  MODULE PROCEDURE bsm_vec_mul_
+END INTERFACE
+
+
+
+CONTAINS
+
+
+
+! ************************************************************************
+! ACCESS PROGRAMS
+! ************************************************************************
+
+
+! ------------------------------------------------------------------------
+!>  \brief Constructor for bandSymMatrixT (non-exception checking)
+!!
+!!  \param self Reference to the banded symmetric matrix to be initialized
+!!  \param hbw  Half bandwidth of the matrix
+!!  \param n    Number of rows/columns in the matrix
+!!
+!!  This routine initializes the bandSymMatrixT object referenced by self.
+!!  The data structures contained in self are allocated and initial values
+!!  are set to zero.
+! ------------------------------------------------------------------------
+SUBROUTINE bsm_init_ (self, hbw,n)
+  TYPE(bandSymMatrixT), INTENT(INOUT) :: self
+  INTEGER, INTENT(IN) :: hbw,n
+  INTEGER :: i,j    !> loop variables
+
+  !> only reallocate if new dimensions do not match existing dimensions
+  IF (bsm_halfBW(self).NE.hbw .OR. bsm_numRows(self).NE.n) THEN
+
+    !> reset object if it was previously initialized
+    CALL bsm_clean(self)
+
+    !> allocate new data structures
+    ALLOCATE(self%dat(hbw,n))
+    ALLOCATE(self%decomp(hbw,n))
+
+  END IF
+
+  !> ensure data is initialized to zero
+  DO j = 1,n
+    DO i = 1,hbw
+      self%dat(i,j) = 0.d0
+      self%decomp(i,j) = 0.d0
+    END DO
+  END DO
+
+  !> reset boolean for decomposition
+  self%is_decomposed = .FALSE.
+
+END SUBROUTINE bsm_init_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Constructor for bandSymMatrixT (exception checking)
+!!
+!!  \param self Reference to the banded symmetric matrix to be initialized
+!!  \param hbw  Half bandwidth of the matrix
+!!  \param n    Number of rows/columns in the matrix
+!!  \param exc  Error code
+!!
+!!  This routine initializes the bandSymMatrixT object referenced by self.
+!!  The data structures contained in self are allocated and initial values
+!!  are set to zero. If the specified dimensions are invalid, a SZE
+!!  exception is returned. If allocation of data structures fails, an
+!!  ALLOC exception is returned.
+! ------------------------------------------------------------------------
+SUBROUTINE bsm_init_exc_ (self, hbw,n, exc)
+  TYPE(bandSymMatrixT), INTENT(INOUT) :: self
+  INTEGER, INTENT(IN) :: hbw,n
+  INTEGER, INTENT(OUT) :: exc
+  INTEGER :: e      !> allocation info code
+  INTEGER :: i,j    !> loop variables
+
+  !> if dimensions are invalid, raise SZE exception
+  IF (hbw.LE.0 .OR. n.LE.0 .OR. hbw.GT.n) THEN
+    exc=SZE
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> only reallocate if new dimensions do not match existing dimensions
+  IF (bsm_halfBW(self).NE.hbw .OR. bsm_numRows(self).NE.n) THEN
+
+    !> reset object if it was previously initialized
+    CALL bsm_clean(self)
+
+    !> allocate new data structure for matrix storage
+    ALLOCATE(self%dat(hbw,n), STAT=e)
+
+    !> if error code is returned, raise ALLOC exception
+    IF (e.NE.0) THEN
+      exc=ALLOC
+      CALL log_printLogMsg(exc,sdr)
+      CALL bsm_clean(self)
+      RETURN
+    ELSE
+      exc=OK
+    END IF
+
+    !> allocate new data structure for decomposition storage
+    ALLOCATE(self%decomp(hbw,n), STAT=e)
+
+    !> if error code is returned, raise ALLOC exception
+    IF (e.NE.0) THEN
+      exc=ALLOC
+      CALL log_printLogMsg(exc,sdr)
+      CALL bsm_clean(self)
+      RETURN
+    ELSE
+      exc=OK
+    END IF
+
+  END IF
+
+  !> ensure data is initialized to zero
+  DO j = 1,n
+    DO i = 1,hbw
+      self%dat(i,j) = 0.d0
+      self%decomp(i,j) = 0.d0
+    END DO
+  END DO
+
+  !> reset boolean for decomposition
+  self%is_decomposed = .FALSE.
+
+END SUBROUTINE bsm_init_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Destructor for bandSymMatrixT
+!!
+!!  \param self Reference to the banded symmetric matrix to be destroyed
+!!
+!!  This routine deallocates existing data structures and resets the
+!!  switch indicating decomposition.
+! ------------------------------------------------------------------------
+SUBROUTINE bsm_clean_ (self)
+  TYPE(bandSymMatrixT), INTENT(INOUT) :: self
+
+  !> deallocate data structures
+  IF (ALLOCATED(self%dat)) DEALLOCATE(self%dat)
+  IF (ALLOCATED(self%decomp)) DEALLOCATE(self%decomp)
+
+  !> reset boolean for decomposition
+  self%is_decomposed = .FALSE.
+
+END SUBROUTINE bsm_clean_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for number of rows
+!!
+!!  \param  self  Reference to the banded symmetric matrix object
+!!  \return n     Number of rows in the matrix
+!!
+!!  This routine determines the number of rows allocated to the matrix
+!!  object. If the matrix is not initialized it returns 0.
+! ------------------------------------------------------------------------
+FUNCTION bsm_num_rows_ (self) RESULT(n)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self
+  INTEGER :: n
+
+  !> if data is initialized, return number of rows
+  IF (ALLOCATED(self%dat)) THEN
+    n = SIZE(self%dat,2)
+  ELSE
+    n = 0   !> if not initialized, there are no rows
+  END IF
+
+END FUNCTION bsm_num_rows_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for half bandwidth
+!!
+!!  \param  self  Reference to the banded symmetric matrix object
+!!  \return hbw   Half bandwidth of the matrix
+!!
+!!  This routine determines the half bandwidth allocated to the matrix
+!!  object. If the matrix is not initialized it returns 0.
+! ------------------------------------------------------------------------
+FUNCTION bsm_half_bw_ (self) RESULT(hbw)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self
+  INTEGER :: hbw
+
+  !> if data is initialized, return half bandwidth
+  IF (ALLOCATED(self%dat)) THEN
+    hbw = SIZE(self%dat,1)
+  ELSE
+    hbw = 0   !> if not initialized, there are no bands
+  END IF
+
+END FUNCTION bsm_half_bw_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for matrix entries (non-exception checking)
+!!
+!!  \param  self  Reference to the banded symmetric matrix object
+!!  \param  i     Row index (in unpacked version)
+!!  \param  j     Column index (in unpacked version)
+!!  \return v     Value at location (i,j) of the matrix
+!!
+!!  This routine determines the value at a particular location in the
+!!  matrix. If the location is outside the band, it returns 0.
+! ------------------------------------------------------------------------
+FUNCTION bsm_get_ (self, i,j) RESULT(v)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self
+  INTEGER, INTENT(IN) :: i,j
+  DOUBLE PRECISION :: v
+  INTEGER :: r,c              !> indices in packed format
+
+  !> get packed indices
+  CALL packed_index(self, i,j, r,c)
+
+  !> if desired location is within the bands, return the value
+  IF (r.GT.0) THEN
+    v = self%dat(r,c)
+  ELSE
+    v = 0.d0    !> if off-band, return 0
+  END IF
+
+END FUNCTION bsm_get_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for matrix entries (exception checking)
+!!
+!!  \param  self  Reference to the banded symmetric matrix object
+!!  \param  i     Row index (in unpacked version)
+!!  \param  j     Column index (in unpacked version)
+!!  \param  exc   Error code
+!!  \return v     Value at location (i,j) of the matrix
+!!
+!!  This routine determines the value at a particular location in the
+!!  matrix. If the location is outside the band, it returns 0. If the
+!!  location is not inside the matrix, it returns a POSIT exception.
+! ------------------------------------------------------------------------
+FUNCTION bsm_get_exc_ (self, i,j, exc) RESULT(v)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self
+  INTEGER, INTENT(IN) :: i,j
+  INTEGER, INTENT(OUT) :: exc
+  DOUBLE PRECISION :: v
+  INTEGER :: n    !> number of rows
+
+  !> get number of rows (for efficiency)
+  n = bsm_numRows(self)
+
+  !> make sure desired indices are inside the matrix
+  IF (i.GT.n .OR. i.LE.0 &
+      .OR. j.GT.n .OR. j.LE.0) THEN
+    exc=POSIT
+    CALL log_printLogMsg(exc,sdr)
+    v = 0.d0
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception getter
+  v = bsm_get(self, i,j)
+
+END FUNCTION bsm_get_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Setter for matrix entries (non-exception checking)
+!!
+!!  \param self Reference to the banded symmetric matrix object
+!!  \param i    Row index (in unpacked version)
+!!  \param j    Column index (in unpacked version)
+!!  \param v    Value to be stored at location (i,j) of the matrix
+!!
+!!  This routine sets the value at a particular location in the
+!!  matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE bsm_set_ (self, i,j,v)
+  TYPE(bandSymMatrixT), INTENT(INOUT) :: self
+  INTEGER, INTENT(IN) :: i,j
+  DOUBLE PRECISION, INTENT(IN) :: v
+  INTEGER :: r,c        !> indices in packed format
+
+  !> get packed indices
+  CALL packed_index(self, i,j, r,c)
+
+  !> set the data value
+  self%dat(r,c) = v
+
+END SUBROUTINE bsm_set_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Setter for matrix entries (exception checking)
+!!
+!!  \param self Reference to the banded symmetric matrix object
+!!  \param i    Row index (in unpacked version)
+!!  \param j    Column index (in unpacked version)
+!!  \param v    Value to be stored at location (i,j) of the matrix
+!!  \param exc  Error code
+!!
+!!  This routine sets the value at a particular location in the
+!!  matrix. If the location is outside the band or if the
+!!  location is not inside the matrix, it returns a POSIT exception.
+! ------------------------------------------------------------------------
+SUBROUTINE bsm_set_exc_ (self, i,j,v, exc)
+  TYPE(bandSymMatrixT), INTENT(INOUT) :: self
+  INTEGER, INTENT(IN) :: i,j
+  DOUBLE PRECISION, INTENT(IN) :: v
+  INTEGER, INTENT(OUT) :: exc
+  INTEGER :: n    !> number of rows
+
+  !> get number of rows (for efficiency)
+  n = bsm_numRows(self)
+
+  !> ensure that desired indices are within the matrix
+  !> and furthermore within the bands
+  IF (i.GT.n .OR. i.LE.0 &
+      .OR. j.GT.n .OR. j.LE.0 &
+      .OR. ABS(j-i).GE.bsm_halfBW(self)) THEN
+    exc=POSIT
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  CALL bsm_set(self, i,j,v)
+
+END SUBROUTINE bsm_set_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Setter for Cholesky decomposition (non-exception checking)
+!!
+!!  \param self   Reference to the banded symmetric matrix object
+!!  \param decomp Array containing the decomposition in packed format
+!!
+!!  This routine sets the Cholesky decomposition of the matrix object.
+!!  The decomposition is typically determined using the DPBTRF subroutine
+!!  from the LAPACK library.
+! ------------------------------------------------------------------------
+SUBROUTINE bsm_set_decomp_ (self,decomp)
+  TYPE(bandSymMatrixT), INTENT(INOUT) :: self
+  DOUBLE PRECISION, INTENT(IN) :: decomp(:,:)
+
+  !> set the decomposition data and flip the boolean
+  self%decomp = decomp
+  self%is_decomposed = .TRUE.
+
+END SUBROUTINE bsm_set_decomp_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Setter for Cholesky decomposition (exception checking)
+!!
+!!  \param self   Reference to the banded symmetric matrix object
+!!  \param decomp Array containing the decomposition in packed format
+!!  \param exc    Error code
+!!
+!!  This routine sets the Cholesky decomposition of the matrix object.
+!!  The decomposition is typically determined using the DPBTRF subroutine
+!!  from the LAPACK library. If the dimensions of the provided
+!!  decomposition array do not match those of the allocated array, it
+!!  returns a DIMEN exception.
+! ------------------------------------------------------------------------
+SUBROUTINE bsm_set_decomp_exc_ (self,decomp, exc)
+  TYPE(bandSymMatrixT), INTENT(INOUT) :: self
+  DOUBLE PRECISION, INTENT(IN) :: decomp(:,:)
+  INTEGER, INTENT(OUT) :: exc
+
+  !> make sure that the dimensions of the proposed decomposition match
+  !> the actual dimensions of the packed storage
+  IF (bsm_halfBW(self).NE.SIZE(decomp,1) &
+        .OR. bsm_numRows(self).NE.SIZE(decomp,2)) THEN
+    exc=DIMEN
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call the non-exception version
+  CALL bsm_setDecomp(self,decomp)
+
+END SUBROUTINE bsm_set_decomp_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for boolean indicating presence of decomposition
+!!
+!!  \param  self          Reference to the banded symmetric matrix object
+!!  \return is_decomposed Boolean indicating whether the decomposition is present
+! ------------------------------------------------------------------------
+FUNCTION bsm_is_decomposed_ (self) RESULT(is_decomposed)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self
+  LOGICAL :: is_decomposed
+
+  is_decomposed = self%is_decomposed
+
+END FUNCTION bsm_is_decomposed_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Add two bandSymMatrixT objects (non-exception checking)
+!!
+!!  \param  self  Reference to the first (left-hand) banded symmetric matrix object
+!!  \param  other Reference to the second (right-hand) banded symmetric matrix object
+!!  \return new   Reference to the resulting banded symmetric matrix object
+!!
+!!  This routine adds two banded symmetric matrix objects. The half
+!!  bandwidth of the result is the maximum of the half bandwidths of the
+!!  input matrix objects.
+! ------------------------------------------------------------------------
+FUNCTION bsm_add_ (self,other) RESULT(new)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self,other
+  TYPE(bandSymMatrixT) :: new
+  INTEGER :: hbw1,hbw2    !> half bandwidths of inputs
+  INTEGER :: hbw,n        !> half bandwidth of output and number of rows
+
+  !> get number of rows and half bandwidths
+  n = bsm_numRows(self)
+  hbw1 = bsm_halfBW(self)
+  hbw2 = bsm_halfBW(other)
+
+  !> new half bandwidth is the maximum of the previous
+  hbw = MAX(hbw1,hbw2)
+
+  !> intialize new bandSymMatrixT
+  CALL bsm_init(new, hbw,n)
+
+  !> add the two matrices (using array ops for efficiency)
+  new%dat(hbw-hbw1+1:hbw,:) = self%dat(:,:)
+  new%dat(hbw-hbw2+1:hbw,:) = new%dat(hbw-hbw2+1:hbw,:) + other%dat(:,:)
+
+END FUNCTION bsm_add_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Add two bandSymMatrixT objects (exception checking)
+!!
+!!  \param  self  Reference to the first (left-hand) banded symmetric matrix object
+!!  \param  other Reference to the second (right-hand) banded symmetric matrix object
+!!  \param  exc   Error code
+!!  \return new   Reference to the resulting banded symmetric matrix object
+!!
+!!  This routine adds two banded symmetric matrix objects. The half
+!!  bandwidth of the result is the maximum of the half bandwidths of the
+!!  input matrix objects. If the dimensions of the two input matrices do
+!!  not match, it returns a DIMEN exception.
+! ------------------------------------------------------------------------
+FUNCTION bsm_add_exc_ (self,other, exc) RESULT(new)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self,other
+  INTEGER, INTENT(OUT) :: exc
+  TYPE(bandSymMatrixT) :: new
+
+  !> make sure that the two matrices have the same number of rows
+  IF (bsm_numRows(self).NE.bsm_numRows(other)) THEN
+    exc=DIMEN
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  new = bsm_add(self,other)
+
+END FUNCTION bsm_add_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Add a matrixT to a bandSymMatrixT in place with mapping (non-exception checking)
+!!
+!!  \param self   Reference to the banded symmetric matrix object
+!!  \param other  Reference to the dense matrix object
+!!  \param ind    Integer mapping indices
+!!
+!!  This routine adds a dense matrix to a banded symmetric matrix in place
+!!  according to an indexed mapping. That is, the element at (i,j) in
+!!  other is added to the element at (ind(i),ind(j)) of self.
+! ------------------------------------------------------------------------
+SUBROUTINE bsm_mapped_add_ (self, other,ind)
+  TYPE(bandSymMatrixT), INTENT(INOUT) :: self
+  TYPE(matrixT), INTENT(IN) :: other
+  INTEGER, INTENT(IN) :: ind(:)
+  INTEGER :: i,j    !> loop variables
+
+  !> Add matrix entries corresponding to mapping indices
+  DO i = 1,SIZE(ind)
+
+    IF (ind(i).EQ.0)  CYCLE
+
+    DO j = i,SIZE(ind)
+
+      IF (ind(j).EQ.0) CYCLE
+
+      CALL bsm_set(self, ind(i),ind(j), &
+                    bsm_get(self,ind(i),ind(j)) &
+                    + dm_get(other,i,j) )
+
+    END DO
+  END DO
+
+END SUBROUTINE bsm_mapped_add_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Add a matrixT to a bandSymMatrixT in place with mapping (exception checking)
+!!
+!!  \param self   Reference to the banded symmetric matrix object
+!!  \param other  Reference to the dense matrix object
+!!  \param ind    Integer mapping indices
+!!  \param exc    Error code
+!!
+!!  This routine adds a dense matrix to a banded symmetric matrix in place
+!!  according to an indexed mapping. That is, the element at (i,j) in
+!!  other is added to the element at (ind(i),ind(j)) of self. If the
+!!  dense matrix is not square or if the its dimensions exceed the half
+!!  bandwidth of the banded symmetric matrix, it returns a DIMEN
+!!  exception. If the number of mapping indices does not match the
+!!  dimension of the dense matrix, it returns a DIMEN exception. If any of
+!!  the index pairs indicate a location outside the matrix bands, it
+!!  returns a POSIT exception. Note that a zero index is valid and
+!!  ind(i)=0 or ind(j)=0 indicate that the value at (i,j) in the dense
+!!  matrix object will not be added to the banded symmetric matrix object.
+! ------------------------------------------------------------------------
+SUBROUTINE bsm_mapped_add_exc_ (self, other,ind, exc)
+  TYPE(bandSymMatrixT), INTENT(INOUT) :: self
+  TYPE(matrixT), INTENT(IN) :: other
+  INTEGER, INTENT(IN) :: ind(:)
+  INTEGER, INTENT(OUT) :: exc
+
+  !> ensure that matrixT is square (it is assumed to be symmetric)
+  !> and that the number of rows in the matrixT matches the number
+  !> of mapping indices
+  !> and that the matrixT is not wider than the half bandwidth
+  !> of the bandSymMatrixT
+  IF (dm_numRows(other).NE.dm_numCols(other) &
+      .OR. dm_numRows(other).NE.SIZE(ind) &
+      .OR. dm_numRows(other).GT.bsm_halfBW(self)) THEN
+    exc=DIMEN
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+
+  !> ensure that there are no mapping indices outside the bandSymMatrixT
+  !> and that all indices are within the bands
+  ELSEIF ( MAXVAL(ind).GT.bsm_numRows(self) .OR. MINVAL(ind).LT.0 &
+            .OR. MAXVAL(ind)-MINVAL(ind,MASK=ind.GT.0) .GE. bsm_halfBW(self) ) THEN
+    exc=POSIT
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  CALL bsm_mappedAdd(self,other,ind)
+
+END SUBROUTINE bsm_mapped_add_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Scalar multiplication (self*k version)
+!!
+!!  \param  self  Reference to the banded symmetric matrix object
+!!  \param  k     Scalar multiplication factor
+!!  \return new   Reference to the resulting banded symmetric matrix object
+!!
+!!  This routine multiplies a banded symmetric matrix by a scalar factor.
+! ------------------------------------------------------------------------
+FUNCTION bsm_scal_mul_mat_scal_ (self,k) RESULT(new)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self
+  DOUBLE PRECISION, INTENT(IN) :: k
+  TYPE(bandSymMatrixT) :: new
+
+  !> initialize output
+  CALL bsm_init(new, bsm_halfBW(self),bsm_numRows(self))
+
+  !> perform multiplication
+  new%dat = k * self%dat
+
+END FUNCTION bsm_scal_mul_mat_scal_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Scalar multiplication (k*self version)
+!!
+!!  \param  self  Reference to the banded symmetric matrix object
+!!  \param  k     Scalar multiplication factor
+!!  \return new   Reference to the resulting banded symmetric matrix object
+!!
+!!  This routine multiplies a banded symmetric matrix by a scalar factor.
+!!  The overload is necessary for the OPERATOR (*) overload. This version
+!!  simply calls the other version for better maintainability.
+! ------------------------------------------------------------------------
+FUNCTION bsm_scal_mul_scal_mat_ (k,self) RESULT(new)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self
+  DOUBLE PRECISION, INTENT(IN) :: k
+  TYPE(bandSymMatrixT) :: new
+
+  !> call self*k version
+  new = bsm_scalMul(self,k)
+
+END FUNCTION bsm_scal_mul_scal_mat_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Post-multiplication of a banded symmetric matrix by a vector (non-exception checking)
+!!
+!!  \param  self  Reference to the banded symmetric matrix object
+!!  \param  other Reference to the vector object
+!!  \return new   Reference to the resulting vector object
+!!
+!!  This routine post-multiplies a banded symmetric matrix by a vector.
+!!  That is, it performs self*other, not other*self.
+! ------------------------------------------------------------------------
+FUNCTION bsm_vec_mul_ (self,other) RESULT(new)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self
+  TYPE(vectorT), INTENT(IN) :: other
+  TYPE(vectorT) :: new
+  INTEGER :: hbw,n    !> half bandwidth and number of rows
+  INTEGER :: i,j      !> loop variables
+
+  !> get size parameters
+  n = vec_length(other)
+  hbw = bsm_halfBW(self)
+
+  !> initialize solution
+  CALL vec_init(new, n)
+
+  !> perform partial dot products of rows in self with the vector
+  !> (ignoring off-band zeros)
+  DO i = 1,n
+    DO j = MAX(i-hbw+1,1), MIN(i+hbw-1,n)
+      CALL vec_set(new, i, vec_get(new,i) + bsm_get(self,i,j)*vec_get(other,j))
+    END DO
+  END DO
+
+END FUNCTION bsm_vec_mul_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Post-multiplication of a banded symmetric matrix by a vector (exception checking)
+!!
+!!  \param  self  Reference to the banded symmetric matrix object
+!!  \param  other Reference to the vector object
+!!  \param  exc   Error code
+!!  \return new   Reference to the resulting vector object
+!!
+!!  This routine post-multiplies a banded symmetric matrix by a vector.
+!!  That is, it performs self*other, not other*self. If the dimensions
+!!  of the two objects are not compatible, it returns a DIMEN exception.
+! ------------------------------------------------------------------------
+FUNCTION bsm_vec_mul_exc_ (self,other, exc) RESULT(new)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self
+  TYPE(vectorT), INTENT(IN) :: other
+  INTEGER, INTENT(OUT) :: exc
+  TYPE(vectorT) :: new
+
+  !> ensure that the number of rows in the bandSymMatrixT and
+  !> the vectorT match
+  IF (bsm_numRows(self).NE.vec_length(other)) THEN
+    exc=DIMEN
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  new = bsm_vecMul(self,other)
+
+END FUNCTION bsm_vec_mul_exc_
+
+
+
+! ************************************************************************
+! LOCAL FUNCTIONS
+! ************************************************************************
+
+
+! ------------------------------------------------------------------------
+!>  \brief Determine the packed indices from the non-packed indices
+!!
+!!  \param self Reference to the banded symmetric matrix object
+!!  \param i    Row index (non-packed version)
+!!  \param j    Column index (non-packed version)
+!!  \param r    Row index (packed version)
+!!  \param c    Column index (packed version)
+!!
+!!  This routine determines the indices of the desired location in the
+!!  matrix in the packed data structure. It accounts for symmetry by
+!!  ensuring that i <= j. After that is checked, the packed indices are
+!!  computed as:  r = hbw - (j-i)
+!!                c = j
+! ------------------------------------------------------------------------
+SUBROUTINE packed_index (self, i,j, r,c)
+  TYPE(bandSymMatrixT), INTENT(IN) :: self
+  INTEGER, INTENT(IN) :: i,j
+  INTEGER, INTENT(OUT) :: r,c
+
+  !> flip indices - symmetry
+  IF (i.GT.j) THEN
+    r=j
+    c=i
+  ELSE
+    r=i
+    c=j
+  END IF
+
+  !> convert to packed index
+  r = bsm_halfBW(self) - (c-r)
+
+END SUBROUTINE packed_index
+
+END MODULE band_sym_matrix_def

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/text-base/band_sym_matrix_test.f90.svn-base

+! ------------------------------------------------------------------------
+!>  \brief Module for testing Banded Symmetric Matrix data type
+! ------------------------------------------------------------------------
+MODULE band_sym_matrix_test
+USE fruit                   !> Unit testing framework
+USE log_message_control     !> Printing log/error messages
+USE log_messages            !> Log/error codes
+USE band_sym_matrix_def     !> Banded Symmetric Matrix data type
+USE dense_matrix_def        !> Dense Matrix data type
+USE vector_def              !> Vector data type
+
+IMPLICIT NONE
+
+
+CONTAINS
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for OK exception message on allocation
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test       Test banded symmetric matrix
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expMsg     Expected log message
+!!  \param actMsg     Actual log message
+!!
+!!  This test initializes a bandSymMatrixT and makes sure that the
+!!  exception message is OK (i.e. allocation did not fail).
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_allocation_MSG
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_allocation_MSG'
+  TYPE(bandSymMatrixT) :: test
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg=OK
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log message file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize matrix and check the exception
+  CALL bsm_init(test, hbw,n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrix
+  CALL log_closeLogFile()
+  CALL bsm_clean(test)
+
+END SUBROUTINE test_bsm_allocation_MSG
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for SZE exception message on allocation
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test banded symmetric matrix
+!!  \param testName Filename for log file (required for exceptions)
+!!  \param hbw      Half bandwidth of the matrix
+!!  \param n        Number of rows/columns in the matrix
+!!  \param expMsg   Expected log message
+!!  \param actMsg   Actual log message
+!!
+!!  This test attempts to initialize a bandSymMatrixT with invalid size
+!!  parameters and verifies that the correct exception is returned.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_allocation_SZE
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_allocation_SZE'
+  TYPE(bandSymMatrixT) :: test
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg=SZE
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> try to initialize with hbw=-1
+  CALL bsm_init(test, -1,n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to initialize with n=-1
+  CALL bsm_init(test, hbw,-1, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrix
+  CALL log_closeLogFile()
+  CALL bsm_clean(test)
+
+END SUBROUTINE test_bsm_allocation_SZE
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for initialization of matrix data to zeros.
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test banded symmetric matrix
+!!  \param hbw      Half bandwidth of the matrix
+!!  \param n        Number of rows/columns in the matrix
+!!  \param expDat   Expected initial contents of matrix data
+!!
+!!  This test initializes a bandSymMatrixT and ensures that the data is
+!!  initialized to zeros.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_allocation_DAT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_allocation_DAT'
+  TYPE(bandSymMatrixT) :: test
+  INTEGER, PARAMETER :: hbw=3,n=10
+  DOUBLE PRECISION, DIMENSION(hbw,n) :: expDat
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up the expected data (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0, &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0, &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0  &
+                    /), SHAPE(expDat) )
+
+  !> initialize matrix and test actual data
+  CALL bsm_init(test, hbw,n)
+  CALL assertEquals(expDat,test%dat, hbw,n)
+  CALL assertEquals(expDat,test%decomp, hbw,n)
+
+  !> deallocate matrix
+  CALL bsm_clean(test)
+
+END SUBROUTINE test_bsm_allocation_DAT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for deallocation
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test       Test banded symmetric matrix
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expBefore  Expected allocation status before deallocation
+!!  \param expAfter   Expected allocation status after deallocation
+!!  \param actBefore  Actual allocation status before deallocation
+!!  \param actAfter   Actual allocation status after deallocation
+!!
+!!  This test initializes a bandSymMatrixT and ensures that it is
+!!  allocated. It then deallocates the object and ensures that it has
+!!  been deallocated.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_deallocation
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_deallocation'
+  TYPE(bandSymMatrixT) :: test
+  INTEGER, PARAMETER :: hbw=3,n=10
+  LOGICAL, PARAMETER :: expBefore = .TRUE.
+  LOGICAL, PARAMETER :: expAfter = .FALSE.
+  LOGICAL :: actBefore, actAfter
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the matrix
+  CALL bsm_init(test, hbw,n)
+
+  !> make sure dat is allocated
+  actBefore = ALLOCATED(test%dat)
+  CALL assertEquals(expBefore,actBefore)
+
+  !> make sure decomp is allocated
+  actBefore = ALLOCATED(test%decomp)
+  CALL assertEquals(expBefore,actBefore)
+
+  !> deallocate the matrix
+  CALL bsm_clean(test)
+
+  !> make sure dat is deallocated
+  actAfter = ALLOCATED(test%dat)
+  CALL assertEquals(expAfter,actAfter)
+
+  !> make sure decomp is deallocated
+  actAfter = ALLOCATED(test%decomp)
+  CALL assertEquals(expAfter,actAfter)
+
+END SUBROUTINE test_bsm_deallocation
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for number of rows when matrix is not allocated
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test       Test banded symmetric matrix
+!!  \param expected   Expected number of rows
+!!  \param actual     Actual number of rows
+!!
+!!  This test makes sure that the number of rows is returned as 0 when
+!!  the matrix is not allocated
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_num_rows_not_allocated
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_num_rows_not_allocated'
+  TYPE(bandSymMatrixT) :: test
+  INTEGER, PARAMETER :: expected = 0
+  INTEGER :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> check number of rows
+  actual = bsm_numRows(test)
+  CALL assertEquals(expected, actual)
+
+END SUBROUTINE test_bsm_num_rows_not_allocated
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for number of rows when matrix is allocated
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test       Test banded symmetric matrix
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expected   Expected number of rows
+!!  \param actual     Actual number of rows
+!!
+!!  This test checks that the number of rows returned is correct when the
+!!  matrix is allocated.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_num_rows_allocated
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_num_rows_allocated'
+  TYPE(bandSymMatrixT) :: test
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expected = n
+  INTEGER :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the matrix
+  CALL bsm_init(test, hbw,n)
+
+  !> check number of rows
+  actual = bsm_numRows(test)
+  CALL assertEquals(expected, actual)
+
+  !> deallocate the matrix
+  CALL bsm_clean(test)
+
+END SUBROUTINE test_bsm_num_rows_allocated
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for half bandwidth when matrix is not allocated
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test       Test banded symmetric matrix
+!!  \param expected   Expected number of rows
+!!  \param actual     Actual number of rows
+!!
+!!  This test makes sure that the half bandwidth is returned as 0 when
+!!  the matrix is not allocated
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_half_bw_not_allocated
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_half_bw_not_allocated'
+  TYPE(bandSymMatrixT) :: test
+  INTEGER, PARAMETER :: expected = 0
+  INTEGER :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> check the half bandwidth
+  actual = bsm_halfBW(test)
+  CALL assertEquals(expected, actual)
+
+END SUBROUTINE test_bsm_half_bw_not_allocated
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for half bandwidth when matrix is allocated
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test       Test banded symmetric matrix
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expected   Expected number of rows
+!!  \param actual     Actual number of rows
+!!
+!!  This test checks that the half bandwidth returned is correct when the
+!!  matrix is allocated.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_half_bw_allocated
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_half_bw_allocated'
+  TYPE(bandSymMatrixT) :: test
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expected = hbw
+  INTEGER :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the matrix
+  CALL bsm_init(test, hbw,n)
+
+  !> check the half bandwidth
+  actual = bsm_halfBW(test)
+  CALL assertEquals(expected, actual)
+
+  !> deallocate the matrix
+  CALL bsm_clean(test)
+
+END SUBROUTINE test_bsm_half_bw_allocated
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for POSIT exception from bsm_get
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test banded symmetric matrix
+!!  \param testName Filename for log file (required for exceptions)
+!!  \param hbw      Half bandwidth of the matrix
+!!  \param n        Number of rows/columns in the matrix
+!!  \param expMsg   Expected log message
+!!  \param actMsg   Actual log message
+!!  \param v        Dummy variable for get function return
+!!
+!!  This test checks that a POSIT exception is returned when the requested
+!!  location is not inside the matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_get_POSIT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_get_POSIT'
+  TYPE(bandSymMatrixT) :: test
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg = POSIT
+  INTEGER :: actMsg
+  DOUBLE PRECISION :: v
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize matrix
+  CALL bsm_init(test, hbw,n)
+
+  !> try to get beyond last row
+  v = bsm_get(test, n+1,n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get before first row
+  v = bsm_get(test, 0,n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get after last column
+  v = bsm_get(test, n,n+1, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get before first column
+  v = bsm_get(test, n,0, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize the log file and deallocate the matrix
+  CALL log_closeLogFile()
+  CALL bsm_clean(test)
+
+END SUBROUTINE test_bsm_get_POSIT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct value return from bsm_get
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test banded symmetric matrix
+!!  \param hbw      Half bandwidth of the matrix
+!!  \param n        Number of rows/columns in the matrix
+!!  \param i        Row index for test location
+!!  \param j        Column index for test location
+!!  \param expVal   Expected data value
+!!  \param actVal   Actual data value
+!!
+!!  This test checks that the correct value is returned from the get
+!!  access program (i.e. that the location is interpreted correctly from
+!!  the unpacked indices).
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_get_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_get_VAL'
+  TYPE(bandSymMatrixT) :: test
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: i=3,j=5
+  DOUBLE PRECISION, PARAMETER :: expVal = 3.d0
+  DOUBLE PRECISION :: actVal
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the matrix
+  CALL bsm_init(test, hbw,n)
+
+  !> set the value manually in the internal data structure
+  test%dat(hbw+i-j,j) = expVal
+
+  !> get the value using the access program
+  actVal = bsm_get(test, i,j)
+  CALL assertEquals(expVal,actVal)
+
+  !> get the value using the access program (testing symmetry)
+  actVal = bsm_get(test, j,i)
+  CALL assertEquals(expVal,actVal)
+
+  !> deallocate the matrix
+  CALL bsm_clean(test)
+
+END SUBROUTINE test_bsm_get_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for POSIT exception from bsm_set
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test banded symmetric matrix
+!!  \param testName Filename for log file (required for exceptions)
+!!  \param hbw      Half bandwidth of the matrix
+!!  \param n        Number of rows/columns in the matrix
+!!  \param expMsg   Expected log message
+!!  \param actMsg   Actual log message
+!!  \param v        Dummy variable for set function input
+!!
+!!  This test checks that a POSIT exception is returned when the requested
+!!  location is not inside the matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_set_POSIT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_set_POSIT'
+  TYPE(bandSymMatrixT) :: test
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg = POSIT
+  DOUBLE PRECISION, PARAMETER :: v=0.d0
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize the matrix
+  CALL bsm_init(test, hbw,n)
+
+  !> try to set after last row
+  CALL bsm_set(test, n+1,n,v, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set before first row
+  CALL bsm_set(test, 0,n,v, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set after last column
+  CALL bsm_set(test, n,n+1,v, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set before first column
+  CALL bsm_set(test, n,0,v, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrix
+  CALL log_closeLogFile()
+  CALL bsm_clean(test)
+
+END SUBROUTINE test_bsm_set_POSIT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct value return after bsm_set
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test banded symmetric matrix
+!!  \param hbw      Half bandwidth of the matrix
+!!  \param n        Number of rows/columns in the matrix
+!!  \param i        Row index for test location
+!!  \param j        Column index for test location
+!!  \param expVal   Expected data value
+!!  \param actVal   Actual data value
+!!  \param reset    For resetting data value before symmetry test
+!!
+!!  This test checks that the correct value is entered by the set
+!!  access program (i.e. that the location is interpreted correctly from
+!!  the unpacked indices).
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_set_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_set_VAL'
+  TYPE(bandSymMatrixT) :: test
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: i=3,j=5
+  DOUBLE PRECISION, PARAMETER :: expVal = 3.d0
+  DOUBLE PRECISION, PARAMETER :: reset = 0.d0
+  DOUBLE PRECISION :: actVal
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the matrix
+  CALL bsm_init(test, hbw,n)
+
+  !> set the value
+  CALL bsm_set(test, i,j,expVal)
+  actVal = bsm_get(test, i,j)
+  CALL assertEquals(expVal,actVal)
+
+  !> set the value (symmetry)
+  CALL bsm_set(test, i,j,reset)
+  CALL bsm_set(test, j,i,expVal)
+  actVal = bsm_get(test, i,j)
+  CALL assertEquals(expVal,actVal)
+
+  !> deallocate matrix
+  CALL bsm_clean(test)
+
+END SUBROUTINE test_bsm_set_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for DIMEN exception on setting decomposed matrix data
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test banded symmetric matrix
+!!  \param testName Filename for log file (required for exceptions)
+!!  \param hbw      Half bandwidth of the matrix
+!!  \param n        Number of rows/columns in the matrix
+!!  \param expMsg   Expected log message
+!!  \param actMsg   Actual log message
+!!  \param decomp   Dummy data that has incorrect dimensions
+!!
+!!  This test checks that a DIMEN exception is returned when the
+!!  decomposed data has incorrect dimensions. That is, not (hbw,n).
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_set_decomp_DIMEN
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_set_decomp_DIMEN'
+  TYPE(bandSymMatrixT) :: test
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg = DIMEN
+  DOUBLE PRECISION, DIMENSION(hbw-1,n) :: decomp
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up dummy data (note: Fortran uses column-major storage)
+  decomp = RESHAPE( (/ &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0, &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0  &
+                    /), SHAPE(decomp) )
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize matrix
+  CALL bsm_init(test, hbw,n)
+
+  !> try to set the decomposed data
+  CALL bsm_setDecomp(test,decomp, actMsg)
+  CALL assertEquals(expMsg, actMsg)
+
+  !> finalize log file and deallocate matrix
+  CALL log_closeLogFile()
+  CALL bsm_clean(test)
+
+END SUBROUTINE test_bsm_set_decomp_DIMEN
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct setting of decomposed data
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test         Test banded symmetric matrix
+!!  \param hbw          Half bandwidth of the matrix
+!!  \param n            Number of rows/columns in the matrix
+!!  \param expDecomp    Expected decomposed data
+!!
+!!  This test checks that the decomposed data is correctly set.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_set_decomp_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_set_decomp_VAL'
+  TYPE(bandSymMatrixT) :: test
+  INTEGER, PARAMETER :: hbw=3,n=10
+  DOUBLE PRECISION, DIMENSION(hbw,n) :: expDecomp
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up decomposed data (note: Fortran uses column-major storage)
+  expDecomp = RESHAPE((/ &
+                      0.d0,0.d0,1.d0,0.d0,2.d0,3.d0,4.d0,5.d0,6.d0,7.d0, &
+                      8.d0,9.d0,0.d0,1.d0,2.d0,3.d0,4.d0,5.d0,6.d0,7.d0, &
+                      8.d0,9.d0,0.d0,1.d0,2.d0,3.d0,4.d0,5.d0,6.d0,7.d0  &
+                      /), SHAPE(expDecomp) )
+
+  !> initialize matrix
+  CALL bsm_init(test, hbw,n)
+
+  !> set the decomposed data
+  CALL bsm_setDecomp(test,expDecomp)
+  CALL assertEquals(expDecomp, test%decomp, hbw,n)
+
+  !> deallocate matrix
+  CALL bsm_clean(test)
+
+END SUBROUTINE test_bsm_set_decomp_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct status of switch for decomposed data
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test       Test banded symmetric matrix
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expBefore  Expected value of switch before setting decomp
+!!  \param expAfter   Expected value of switch after setting decomp
+!!  \param expClean   Expected value of switch after deallocating matrix
+!!  \param actBefore  Actual value of switch before setting decomp
+!!  \param actAfter   Actual value of switch after setting decomp
+!!  \param actClean   Actual value of switch after deallocating matrix
+!!  \param decomp     Dummy decomposed data
+!!
+!!  This test checks that the presence of decomposed data is correctly
+!!  identified.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_is_decomposed
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_is_decomposed'
+  TYPE(bandSymMatrixT) :: test
+  INTEGER, PARAMETER :: hbw=3,n=10
+  LOGICAL, PARAMETER :: expBefore = .FALSE.
+  LOGICAL, PARAMETER :: expAfter = .TRUE.
+  LOGICAL, PARAMETER :: expClean = .FALSE.
+  DOUBLE PRECISION, DIMENSION(hbw,n) :: decomp
+  LOGICAL :: actBefore, actAfter, actClean
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize matrix and check initial value of switch
+  CALL bsm_init(test, hbw,n)
+  actBefore = bsm_isDecomposed(test)
+  CALL assertEquals(expBefore,actBefore)
+
+  !> set up decomposed data (note: Fortran uses column-major storage)
+  decomp = RESHAPE( (/ &
+                    0.d0,0.d0,1.d0,0.d0,2.d0,3.d0,4.d0,5.d0,6.d0,7.d0, &
+                    8.d0,9.d0,0.d0,1.d0,2.d0,3.d0,4.d0,5.d0,6.d0,7.d0, &
+                    8.d0,9.d0,0.d0,1.d0,2.d0,3.d0,4.d0,5.d0,6.d0,7.d0  &
+                    /), SHAPE(decomp) )
+
+  !> set decomposed data
+  CALL bsm_setDecomp(test,decomp)
+  actAfter = bsm_isDecomposed(test)
+  CALL assertEquals(expAfter,actAfter)
+
+  !> deallocate matrix and check final value of decomposed switch
+  CALL bsm_clean(test)
+  actClean = bsm_isDecomposed(test)
+  CALL assertEquals(expClean,actClean)
+
+END SUBROUTINE test_bsm_is_decomposed
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for DIMEN exception when adding matrices
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      First input banded symmetric matrix
+!!  \param test2      Second input banded symmetric matrix
+!!  \param test3      Output banded symmetric matrix
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expMsg     Expected exception message
+!!  \param actMsg     Actual exception message
+!!
+!!  This test checks that a DIMEN exception is raised when matrices with
+!!  incompatible dimensions are added.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_add_DIMEN
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_add_DIMEN'
+  TYPE(bandSymMatrixT) :: test1, test2, test3
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg=DIMEN
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize input matrices (data contents not important)
+  CALL bsm_init(test1, hbw,n)
+  CALL bsm_init(test2, hbw,n+1) !> note different dimension
+
+  !> attempt to add the matrices
+  test3 = bsm_add(test1,test2, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrices
+  CALL log_closeLogFile()
+  CALL bsm_clean(test1)
+  CALL bsm_clean(test2)
+  CALL bsm_clean(test3)
+
+END SUBROUTINE test_bsm_add_DIMEN
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct addition of matrices
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      First input banded symmetric matrix
+!!  \param test2      Second input banded symmetric matrix
+!!  \param test3      Output banded symmetric matrix
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expDat     Expected result of addition
+!!
+!!  This test checks that the correct result is obtained when adding two
+!!  banded symmetric matrices.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_add_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_add_VAL'
+  TYPE(bandSymMatrixT) :: test1, test2, test3
+  INTEGER, PARAMETER :: hbw=3,n=10
+  DOUBLE PRECISION, DIMENSION(hbw,n) :: expDat
+  INTEGER :: i,j  !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize expected result (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    0.d0,0.d0,3.d0,0.d0,3.d0,3.d0,1.d0,3.d0,3.d0,1.d0, &
+                    3.d0,3.d0,1.d0,3.d0,3.d0,1.d0,3.d0,3.d0,1.d0,3.d0, &
+                    3.d0,1.d0,3.d0,3.d0,1.d0,3.d0,3.d0,1.d0,3.d0,3.d0  &
+                    /), SHAPE(expDat) )
+
+  !> initialize first input matrix
+  !! (loop sets matrix on bands to all ones)
+  CALL bsm_init(test1, hbw,n)
+  DO i = 1,n
+    DO j = i,MIN(i+hbw-1,n)
+      CALL bsm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> initialize second input matrix (note different bandwidth)
+  !! (loop sets matrix on bands to all twos)
+  CALL bsm_init(test2, hbw-1,n)
+  DO i = 1,n
+    DO j = i,MIN(i+(hbw-1)-1,n)
+     CALL bsm_set(test2, i,j, 2.d0)
+    END DO
+  END DO
+
+  !> add the matrices
+  test3 = bsm_add(test1,test2)
+  CALL assertEquals(expDat,test3%dat, hbw,n)
+
+  !> deallocate matrices
+  CALL bsm_clean(test1)
+  CALL bsm_clean(test2)
+  CALL bsm_clean(test3)
+
+END SUBROUTINE test_bsm_add_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct addition of matrices ( with OPERATOR (+) )
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      First input banded symmetric matrix
+!!  \param test2      Second input banded symmetric matrix
+!!  \param test3      Output banded symmetric matrix
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expDat     Expected result of addition
+!!
+!!  This test checks that the correct result is obtained when adding two
+!!  banded symmetric matrices using the (+) operator.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_add_OP
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_add_OP'
+  TYPE(bandSymMatrixT) :: test1, test2, test3
+  INTEGER, PARAMETER :: hbw=3,n=10
+  DOUBLE PRECISION, DIMENSION(hbw,n) :: expDat
+  INTEGER :: i,j  !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize expected result (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    0.d0,0.d0,3.d0,0.d0,3.d0,3.d0,1.d0,3.d0,3.d0,1.d0, &
+                    3.d0,3.d0,1.d0,3.d0,3.d0,1.d0,3.d0,3.d0,1.d0,3.d0, &
+                    3.d0,1.d0,3.d0,3.d0,1.d0,3.d0,3.d0,1.d0,3.d0,3.d0  &
+                    /), SHAPE(expDat) )
+
+  !> initialize first input matrix
+  !! (loop sets matrix on bands to all ones)
+  CALL bsm_init(test1, hbw,n)
+  DO i = 1,n
+    DO j = i,MIN(i+hbw-1,n)
+      CALL bsm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> initialize second input matrix (note different bandwidth)
+  !! (loop sets matrix on bands to all twos)
+  CALL bsm_init(test2, hbw-1,n)
+  DO i = 1,n
+    DO j = i,MIN(i+(hbw-1)-1,n)
+     CALL bsm_set(test2, i,j, 2.d0)
+    END DO
+  END DO
+
+  !> add the matrices
+  test3 = test1 + test2
+  CALL assertEquals(expDat,test3%dat, hbw,n)
+
+  !> deallocate matrices
+  CALL bsm_clean(test1)
+  CALL bsm_clean(test2)
+  CALL bsm_clean(test3)
+
+END SUBROUTINE test_bsm_add_OP
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for DIMEN exception when performing mapped add
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Test banded symmetric matrix
+!!  \param test2      Test input dense matrix
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expMsg     Expected exception message
+!!  \param actMsg     Actual exception message
+!!  \param ind        Mapping indices (size does not match test2)
+!!
+!!  This test checks that a DIMEN exception is raised when mapped addition
+!!  is attempted with the dimension of the index array not matching the
+!!  dimensions of the dense matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_mapped_add_DIMEN1
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_mapped_add_DIMEN1'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(matrixT) :: test2
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg=DIMEN
+  INTEGER, PARAMETER, DIMENSION(hbw-1) :: ind = (/ 1,2 /)
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize matrices (data content not important)
+  CALL bsm_init(test1, hbw,n)
+  CALL dm_init(test2, hbw,hbw)
+
+  !> try mapped add (note hbw=3, SIZE(ind)=2)
+  CALL bsm_mappedAdd(test1, test2,ind, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrices
+  CALL log_closeLogFile()
+  CALL bsm_clean(test1)
+  CALL dm_clean(test2)
+
+END SUBROUTINE test_bsm_mapped_add_DIMEN1
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for DIMEN exception when performing mapped add
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Test banded symmetric matrix
+!!  \param test2      Test input dense matrix (not square)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expMsg     Expected exception message
+!!  \param actMsg     Actual exception message
+!!  \param ind        Mapping indices
+!!
+!!  This test checks that a DIMEN exception is raised when mapped addition
+!!  is attempted with a rectangular input dense matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_mapped_add_DIMEN2
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_mapped_add_DIMEN2'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(matrixT) :: test2
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg=DIMEN
+  INTEGER, PARAMETER, DIMENSION(hbw) :: ind = (/ 1,2,3 /)
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize matrices (note test2 is not square)
+  CALL bsm_init(test1, hbw,n)
+  CALL dm_init(test2, hbw+1,hbw)
+
+  !> attempt mapped add
+  CALL bsm_mappedAdd(test1, test2,ind, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrices
+  CALL log_closeLogFile()
+  CALL bsm_clean(test1)
+  CALL dm_clean(test2)
+
+END SUBROUTINE test_bsm_mapped_add_DIMEN2
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for DIMEN exception when performing mapped add
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Test banded symmetric matrix
+!!  \param test2      Test input dense matrix (larger than hbw)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expMsg     Expected exception message
+!!  \param actMsg     Actual exception message
+!!  \param ind        Mapping indices (larger than hbw)
+!!
+!!  This test checks that a DIMEN exception is raised when mapped addition
+!!  is attempted with an input dense matrix whose size exceeds the half
+!!  bandwidth of the banded symmetric matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_mapped_add_DIMEN3
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_mapped_add_DIMEN3'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(matrixT) :: test2
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg=DIMEN
+  INTEGER, PARAMETER, DIMENSION(hbw+1) :: ind = (/ 1,2,3,4 /)
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize matrices (note that test2 is wider than hbw)
+  CALL bsm_init(test1, hbw,n)
+  CALL dm_init(test2, hbw+1,hbw+1)
+
+  !> attempt mapped add
+  CALL bsm_mappedAdd(test1, test2,ind, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrices
+  CALL log_closeLogFile()
+  CALL bsm_clean(test1)
+  CALL dm_clean(test2)
+
+END SUBROUTINE test_bsm_mapped_add_DIMEN3
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for POSIT exception when performing mapped add
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Test banded symmetric matrix
+!!  \param test2      Test input dense matrix
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expMsg     Expected exception message
+!!  \param actMsg     Actual exception message
+!!  \param ind1       Mapping indices (with index below minimum)
+!!  \param ind2       Mapping indices (with index above maximum)
+!!
+!!  This test checks that a POSIT exception is raised when mapped addition
+!!  is attempted with mapping indices that are outside the matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_mapped_add_POSIT1
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_mapped_add_POSIT1'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(matrixT) :: test2
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg=POSIT
+  INTEGER, PARAMETER, DIMENSION(hbw) :: ind1 = (/ -1,0,1 /)
+  INTEGER, PARAMETER, DIMENSION(hbw) :: ind2 = (/ 9,10,11 /)
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize matrices (data not important)
+  CALL bsm_init(test1, hbw,n)
+  CALL dm_init(test2, hbw,hbw)
+
+  !> attempt first mapped add
+  CALL bsm_mappedAdd(test1, test2,ind1, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> attempt second mapped add
+  CALL bsm_mappedAdd(test1, test2,ind2, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrices
+  CALL log_closeLogFile()
+  CALL bsm_clean(test1)
+  CALL dm_clean(test2)
+
+END SUBROUTINE test_bsm_mapped_add_POSIT1
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for POSIT exception when performing mapped add
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Test banded symmetric matrix
+!!  \param test2      Test input dense matrix
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expMsg     Expected exception message
+!!  \param actMsg     Actual exception message
+!!  \param ind        Mapping indices (range of indices >= half bandwidth)
+!!
+!!  This test checks that a POSIT exception is raised when mapped addition
+!!  is attempted with mapping indices that have a difference equal to the
+!!  half bandwidth (i.e. trying to set a value outside the band).
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_mapped_add_POSIT2
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_mapped_add_POSIT2'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(matrixT) :: test2
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg=POSIT
+  INTEGER, PARAMETER, DIMENSION(hbw) :: ind = (/ 4,6,9 /)
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize matrices (data content not important)
+  CALL bsm_init(test1, hbw,n)
+  CALL dm_init(test2, hbw,hbw)
+
+  !> attempt mapped add
+  CALL bsm_mappedAdd(test1, test2,ind, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrices
+  CALL log_closeLogFile()
+  CALL bsm_clean(test1)
+  CALL dm_clean(test2)
+
+END SUBROUTINE test_bsm_mapped_add_POSIT2
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct mapped addition (single operation)
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Test banded symmetric matrix
+!!  \param test2      Test input dense matrix
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expDat     Expected data contents
+!!  \param ind        Mapping indices
+!!
+!!  This test checks that mapped addition is performed correctly. The
+!!  operation is performed once on a blank matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_mapped_add_VAL1
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_mapped_add_VAL1'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(matrixT) :: test2
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER, DIMENSION(2) :: ind = (/ 1,3 /)
+  DOUBLE PRECISION, DIMENSION(hbw,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set expected result of mapped addition (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    0.d0,0.d0,1.d0,0.d0,0.d0,0.d0,2.d0,0.d0,3.d0,0.d0, &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0, &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0  &
+                    /), SHAPE(expDat) )
+
+  !> initialize target matrix
+  CALL bsm_init(test1, hbw,n)
+
+  !> initialize input matrix
+  !! (loop sets the data contents)
+  CALL dm_init(test2, 2,2)
+  DO j = 1,2
+    DO i = 1,2
+      CALL dm_set(test2, i,j, ( (DBLE(j)-1.d0) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> perform mapped addition
+  CALL bsm_mappedAdd(test1, test2,ind)
+  CALL assertEquals(expDat,test1%dat, hbw,n)
+
+  !> deallocate matrices
+  CALL bsm_clean(test1)
+  CALL dm_clean(test2)
+
+END SUBROUTINE test_bsm_mapped_add_VAL1
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct mapped addition (zero index)
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Test banded symmetric matrix
+!!  \param test2      Test input dense matrix
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expDat     Expected data contents
+!!  \param ind        Mapping indices (with zero index indicating no-op)
+!!
+!!  This test checks that mapped addition is performed correctly. The
+!!  mapping index contains a zero, which should tell the mapped adder to
+!!  skip that location in the input matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_mapped_add_VAL2
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_mapped_add_VAL2'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(matrixT) :: test2
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER, DIMENSION(3) :: ind = (/ 1,0,3 /)
+  DOUBLE PRECISION, DIMENSION(hbw,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    0.d0,0.d0,1.d0,0.d0,0.d0,0.d0,3.d0,0.d0,5.d0,0.d0, &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0, &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0  &
+                    /), SHAPE(expDat) )
+
+  !> initialize target matrix
+  CALL bsm_init(test1, hbw,n)
+
+  !> initialize input matrix
+  !! (loop sets the input data)
+  CALL dm_init(test2, 3,3)
+  DO j = 1,3
+    DO i = 1,3
+      CALL dm_set(test2, i,j, ( (DBLE(j)-1.d0) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> perform mapped addition
+  CALL bsm_mappedAdd(test1, test2,ind)
+  CALL assertEquals(expDat,test1%dat, hbw,n)
+
+  !> deallocate matrices
+  CALL bsm_clean(test1)
+  CALL dm_clean(test2)
+
+END SUBROUTINE test_bsm_mapped_add_VAL2
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct mapped addition (multiple operations)
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Test banded symmetric matrix
+!!  \param test2      First test input dense matrix
+!!  \param test3      Second test input dense matrix
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expDat     Expected data contents
+!!  \param ind2       Mapping indices for test2
+!!  \param ind3       Mapping indices for test3
+!!
+!!  This test checks that mapped addition is performed correctly. The
+!!  operation is performed twice to ensure summation of overlapping
+!!  values. That is, input is added to existing data in the target matrix
+!!  rather than overwriting.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_mapped_add_VAL3
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_mapped_add_VAL3'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(matrixT) :: test2,test3
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER, DIMENSION(2) :: ind2 = (/ 1,3 /)
+  INTEGER, PARAMETER, DIMENSION(3) :: ind3 = (/ 1,0,3 /)
+  DOUBLE PRECISION, DIMENSION(hbw,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    0.d0,0.d0,2.d0,0.d0,0.d0,0.d0,5.d0,0.d0,8.d0,0.d0, &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0, &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0  &
+                    /), SHAPE(expDat) )
+
+  !> initialize target matrix
+  CALL bsm_init(test1, hbw,n)
+
+  !> initialize first input matrix
+  !! (loop sets up the data)
+  CALL dm_init(test2, 2,2)
+  DO j = 1,2
+    DO i = 1,2
+      CALL dm_set(test2, i,j, ( (DBLE(j)-1.d0) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> initialize second input matrix
+  !! (loop sets up the data)
+  CALL dm_init(test3, 3,3)
+  DO j = 1,3
+    DO i = 1,3
+      CALL dm_set(test3, i,j, ( (DBLE(j)-1.d0) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> perform the mapped addition for both inputs
+  CALL bsm_mappedAdd(test1, test2,ind2)
+  CALL bsm_mappedAdd(test1, test3,ind3)
+  CALL assertEquals(expDat,test1%dat, hbw,n)
+
+  !> deallocate matrices
+  CALL bsm_clean(test1)
+  CALL dm_clean(test2)
+  CALL dm_clean(test3)
+
+END SUBROUTINE test_bsm_mapped_add_VAL3
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for scalar multiplication
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1    Input banded symmetric matrix
+!!  \param test2    Output banded symmetric matrix
+!!  \param k        Scalar coefficient
+!!  \param hbw      Half bandwidth of the matrix
+!!  \param n        Number of rows/columns in the matrix
+!!  \param expDat   Expected data contents
+!!
+!!  This test checks that scalar multiplication is performed correctly.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_scal_mul_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_scal_mul_VAL'
+  TYPE(bandSymMatrixT) :: test1, test2
+  DOUBLE PRECISION, PARAMETER :: k=3.d0
+  INTEGER, PARAMETER :: hbw=3,n=10
+  DOUBLE PRECISION, DIMENSION(hbw,n) :: expDat
+  INTEGER :: i,j    !> loop variable
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    0.d0,0.d0,3.d0,0.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0, &
+                    3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0, &
+                    3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0  &
+                    /), SHAPE(expDat) )
+
+  !> initialize matrix
+  !! (loop sets band to ones)
+  CALL bsm_init(test1, hbw,n)
+  DO i = 1,n
+    DO j = i,MIN(i+hbw-1,n)
+      CALL bsm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> multiply by scalar
+  test2 = bsm_scalMul(test1,k)
+  CALL assertEquals(expDat,test2%dat, hbw,n)
+
+  !> deallocate matrices
+  CALL bsm_clean(test1)
+  CALL bsm_clean(test2)
+
+END SUBROUTINE test_bsm_scal_mul_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for scalar multiplication (by zero)
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1    Input banded symmetric matrix
+!!  \param test2    Output banded symmetric matrix
+!!  \param k        Scalar coefficient (zero)
+!!  \param hbw      Half bandwidth of the matrix
+!!  \param n        Number of rows/columns in the matrix
+!!  \param expDat   Expected data contents
+!!
+!!  This test checks that scalar multiplication results in a zero matrix
+!!  when the scalar factor is zero.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_scal_mul_ZERO
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_scal_mul_ZERO'
+  TYPE(bandSymMatrixT) :: test1, test2
+  DOUBLE PRECISION, PARAMETER :: k=0.d0
+  INTEGER, PARAMETER :: hbw=3,n=10
+  DOUBLE PRECISION, DIMENSION(hbw,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected data (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0, &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0, &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0  &
+                    /), SHAPE(expDat) )
+
+  !> initialize matrix
+  !! (loop sets band to ones)
+  CALL bsm_init(test1, hbw,n)
+  DO i = 1,n
+    DO j = i,MIN(i+hbw-1,n)
+      CALL bsm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> multiply by zero
+  test2 = bsm_scalMul(test1,k)
+  CALL assertEquals(expDat,test2%dat, hbw,n)
+
+  !> deallocate matrices
+  CALL bsm_clean(test1)
+  CALL bsm_clean(test2)
+
+END SUBROUTINE test_bsm_scal_mul_ZERO
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for scalar multiplication ( using OPERATOR (*) )
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1    Input banded symmetric matrix
+!!  \param test2    Output banded symmetric matrix
+!!  \param k        Scalar coefficient
+!!  \param hbw      Half bandwidth of the matrix
+!!  \param n        Number of rows/columns in the matrix
+!!  \param expDat   Expected data contents
+!!
+!!  This test checks that scalar multiplication is performed correctly
+!!  when using the overloaded (*) operator.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_scal_mul_OP
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_scal_mul_OP'
+  TYPE(bandSymMatrixT) :: test1, test2
+  DOUBLE PRECISION, PARAMETER :: k=3.d0
+  INTEGER, PARAMETER :: hbw=3,n=10
+  DOUBLE PRECISION, DIMENSION(hbw,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    0.d0,0.d0,3.d0,0.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0, &
+                    3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0, &
+                    3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0,3.d0  &
+                    /), SHAPE(expDat) )
+
+  !> initialize input matrix
+  !! (loop sets bands to ones)
+  CALL bsm_init(test1, hbw,n)
+  DO i = 1,n
+    DO j = i,MIN(i+hbw-1,n)
+      CALL bsm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> perform scalar multiplication
+  test2 = k * test1
+  CALL assertEquals(expDat,test2%dat, hbw,n)
+
+  !> clear output matrix
+  CALL bsm_clean(test2)
+
+  !> perform scalar multiplication (switch order of inputs)
+  test2 = test1 * k
+  CALL assertEquals(expDat,test2%dat, hbw,n)
+
+  !> deallocate matrices
+  CALL bsm_clean(test1)
+  CALL bsm_clean(test2)
+
+END SUBROUTINE test_bsm_scal_mul_OP
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for DIMEN exception when performing matrix-vector multiplication
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Input banded symmetric matrix
+!!  \param test2      Input vector
+!!  \param test3      Dummy output vector
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expMsg     Expected exception message
+!!  \param actMsg     Actual exception message
+!!
+!!  This test checks that a DIMEN exception is raised when matrix-vector
+!!  multiplication is attempted with a vector whose length does not equal
+!!  the number of rows in the matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_vec_mul_DIMEN
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_vec_mul_DIMEN'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(vectorT) :: test2, test3
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: hbw=3,n=10
+  INTEGER, PARAMETER :: expMsg=DIMEN
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> set up matrix and vector (data contents not important)
+  CALL bsm_init(test1, hbw,n)
+  CALL vec_init(test2, n+1)   !> note different size
+
+  !> attempt matrix-vector multiplication
+  test3 = bsm_vecMul(test1,test2, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate objects
+  CALL log_closeLogFile()
+  CALL bsm_clean(test1)
+  CALL vec_clean(test2)
+  CALL vec_clean(test3)
+
+END SUBROUTINE test_bsm_vec_mul_DIMEN
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct matrix-vector multiplication
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Input banded symmetric matrix
+!!  \param test2      Input vector
+!!  \param test3      Output vector
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-vector multiplication is performed
+!!  correctly in the general case (input matrix and vector contain data,
+!!  but it is not of a special form).
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_vec_mul_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_vec_mul_VAL'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(vectorT) :: test2, test3
+  INTEGER, PARAMETER :: hbw=3,n=10
+  DOUBLE PRECISION, DIMENSION(n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set expected result (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    14.d0,  40.d0,  85.d0,  150.d0, 235.d0, &
+                    340.d0, 465.d0, 610.d0, 566.d0, 488.d0  &
+                    /), SHAPE(expDat) )
+
+  !> initialize input matrix
+  !! (loop sets up data)
+  CALL bsm_init(test1, hbw,n)
+  DO i = 1,n
+    DO j = i,MIN(i+hbw-1,n)
+      CALL bsm_set(test1, i,j, ( (DBLE(j)-1.d0) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> initialize input vector
+  !! (loop sets up data)
+  CALL vec_init(test2, n)
+  DO i = 1,n
+    CALL vec_set(test2, i, DBLE(i))
+  END DO
+
+  !> perform multiplication
+  test3 = bsm_vecMul(test1,test2)
+  CALL assertEquals(expDat,test3%dat, n)
+
+  !> deallocate objects
+  CALL bsm_clean(test1)
+  CALL vec_clean(test2)
+  CALL vec_clean(test3)
+
+END SUBROUTINE test_bsm_vec_mul_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for matrix-vector multiplication (by the zero vector)
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Input banded symmetric matrix
+!!  \param test2      Input vector (zeros)
+!!  \param test3      Output vector
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-vector multiplication is performed
+!!  correctly when the input vector is all zeros.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_vec_mul_ZERO
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_vec_mul_ZERO'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(vectorT) :: test2, test3
+  INTEGER, PARAMETER :: hbw=3,n=10
+  DOUBLE PRECISION, DIMENSION(n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0,0.d0 &
+                    /), SHAPE(expDat) )
+
+  !> initialize input matrix
+  !! (loop sets up data)
+  CALL bsm_init(test1, hbw,n)
+  DO i = 1,n
+    DO j = i,MIN(i+hbw-1,n)
+      CALL bsm_set(test1, i,j, ( (DBLE(j)-1.d0) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> initialize input vector (all zeros)
+  CALL vec_init(test2, n)
+
+  !> perform multiplication
+  test3 = bsm_vecMul(test1,test2)
+  CALL assertEquals(expDat,test3%dat, n)
+
+  !> deallocate objects
+  CALL bsm_clean(test1)
+  CALL vec_clean(test2)
+  CALL vec_clean(test3)
+
+END SUBROUTINE test_bsm_vec_mul_ZERO
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for matrix-vector multiplication (by the identity matrix)
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Input banded symmetric matrix (identity)
+!!  \param test2      Input vector
+!!  \param test3      Output vector
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-vector multiplication is performed
+!!  correctly when the input matrix is the identity for multiplication.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_vec_mul_IDENT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_vec_mul_IDENT'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(vectorT) :: test2, test3
+  INTEGER, PARAMETER :: hbw=3,n=10
+  DOUBLE PRECISION, DIMENSION(n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize input matrix
+  !! (loop sets up identity matrix)
+  CALL bsm_init(test1, hbw,n)
+  DO i = 1,n
+    CALL bsm_set(test1, i,i, 1.d0)
+  END DO
+
+  !> initialize input vector
+  !! (loop sets data)
+  CALL vec_init(test2, n)
+  DO i = 1,n
+    CALL vec_set(test2, i, DBLE(i))
+  END DO
+
+  !> expected result is the same as the initial vector
+  expDat = test2%dat
+
+  !> perform multiplication
+  test3 = bsm_vecMul(test1,test2)
+  CALL assertEquals(expDat,test3%dat, n)
+
+  !> deallocate objects
+  CALL bsm_clean(test1)
+  CALL vec_clean(test2)
+  CALL vec_clean(test3)
+
+END SUBROUTINE test_bsm_vec_mul_IDENT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct matrix-vector multiplication ( using OPERATOR (*) )
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1      Input banded symmetric matrix
+!!  \param test2      Input vector
+!!  \param test3      Output vector
+!!  \param hbw        Half bandwidth of the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-vector multiplication is performed
+!!  correctly in the general case (input matrix and vector contain data,
+!!  but it is not of a special form). The overloaded operator (*) is used.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bsm_vec_mul_OP
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bsm_vec_mul_OP'
+  TYPE(bandSymMatrixT) :: test1
+  TYPE(vectorT) :: test2, test3
+  INTEGER, PARAMETER :: hbw=3,n=10
+  DOUBLE PRECISION, DIMENSION(n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result (note: Fortran uses column-major storage)
+  expDat = RESHAPE( (/ &
+                    14.d0,  40.d0,  85.d0,  150.d0, 235.d0, &
+                    340.d0, 465.d0, 610.d0, 566.d0, 488.d0  &
+                    /), SHAPE(expDat) )
+
+  !> initialize input matrix
+  !! (loop sets data)
+  CALL bsm_init(test1, hbw,n)
+  DO i = 1,n
+    DO j = i,MIN(i+hbw-1,n)
+      CALL bsm_set(test1, i,j, ( (DBLE(j)-1.d0) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> initialize input vector
+  !! (loop sets data)
+  CALL vec_init(test2, n)
+  DO i = 1,n
+    CALL vec_set(test2, i, DBLE(i))
+  END DO
+
+  !> perform multiplication
+  test3 = test1 * test2
+  CALL assertEquals(expDat,test3%dat, n)
+
+  !> deallocate objects
+  CALL bsm_clean(test1)
+  CALL vec_clean(test2)
+  CALL vec_clean(test3)
+
+END SUBROUTINE test_bsm_vec_mul_OP
+
+
+
+END MODULE band_sym_matrix_test

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/text-base/boundary_data.f90.svn-base

+! ------------------------------------------------------------------------
+!>  \brief Module for Boundary Data
+! ------------------------------------------------------------------------
+MODULE boundary_data
+USE system_constants      !> Global system constants
+USE log_message_control   !> Printing log/error messages
+USE log_messages          !> Log/error message and sender codes
+USE field_data            !> Field Data module
+
+IMPLICIT NONE
+PRIVATE
+
+
+
+! ************************************************************************
+! EXPORTS
+! ************************************************************************
+
+
+!> Exported data types
+PUBLIC :: surfLoadT
+
+
+!> Exported interfaces
+PUBLIC :: bnd_init, bnd_clean, &
+          bnd_numBoundElem, &
+          bnd_getConnect, bnd_setConnect, &
+          bnd_lenBoundElem, &
+          bnd_getTrac, bnd_setTrac
+
+
+
+! ************************************************************************
+! LOCAL CONSTANTS
+! ************************************************************************
+
+!> sender code for this module
+INTEGER, PARAMETER :: sdr = BNDDAT
+
+
+
+! ************************************************************************
+! DATA TYPES
+! ************************************************************************
+
+
+! ------------------------------------------------------------------------
+!>  \brief Surface load
+!!
+!!  \param sig_nt   Shear stress
+!!  \param sig_nn   Normal stress
+! ------------------------------------------------------------------------
+TYPE surfLoadT
+  DOUBLE PRECISION :: sig_nt, sig_nn
+END TYPE surfLoadT
+
+
+! ------------------------------------------------------------------------
+!>  \brief Surface traction
+!!
+!!  \param loads   Set of node loads
+! ------------------------------------------------------------------------
+TYPE tracT
+  TYPE(surfLoadT), DIMENSION(NNODELB) :: loads
+END TYPE tracT
+
+
+! ------------------------------------------------------------------------
+!>  \brief Connectivity
+!!
+!!  \param nd   Sequence of node indices
+! ------------------------------------------------------------------------
+TYPE boundConnectT
+  INTEGER, DIMENSION(NNODELB) :: nd
+END TYPE boundConnectT
+
+
+! ------------------------------------------------------------------------
+!>  \brief Traction Element
+!!
+!!  \param num      Element number
+!!  \param ico      Connectivity of local nodes
+!!  \param trac     Surface traction of element
+! ------------------------------------------------------------------------
+TYPE boundElementT
+  INTEGER             :: num
+  TYPE(boundConnectT) :: ico
+  TYPE(tracT)         :: trac
+END TYPE boundElementT
+
+
+
+! ************************************************************************
+! STATE VARIABLES
+! ************************************************************************
+
+TYPE(boundElementT), ALLOCATABLE :: boundElements(:)  !> set of traction elements
+
+
+
+! ************************************************************************
+! INTERFACES
+! ************************************************************************
+
+
+!> \brief Interface to initializer for boundary element data
+INTERFACE bnd_init
+  MODULE PROCEDURE bnd_init_
+  MODULE PROCEDURE bnd_init_exc_
+END INTERFACE bnd_init
+
+!> \brief Interface to destructor for boundary element data
+INTERFACE bnd_clean
+  MODULE PROCEDURE bnd_clean_
+END INTERFACE bnd_clean
+
+!> \brief Interface to getter for number of traction elements
+INTERFACE bnd_numBoundElem
+  MODULE PROCEDURE bnd_num_bound_elem_
+END INTERFACE bnd_numBoundElem
+
+!> \brief Interface to getter for traction element connectivity
+INTERFACE bnd_getConnect
+  MODULE PROCEDURE bnd_get_connect_
+  MODULE PROCEDURE bnd_get_connect_exc_
+END INTERFACE bnd_getConnect
+
+!> \brief Interface to setter for traction element connectivity
+INTERFACE bnd_setConnect
+  MODULE PROCEDURE bnd_set_connect_
+  MODULE PROCEDURE bnd_set_connect_exc_
+END INTERFACE bnd_setConnect
+
+!> \brief Interface to getter for length of traction element
+INTERFACE bnd_lenBoundElem
+  MODULE PROCEDURE bnd_len_bound_elem_
+  MODULE PROCEDURE bnd_len_bound_elem_exc_
+END INTERFACE bnd_lenBoundElem
+
+!> \brief Interface to getter for tractions
+INTERFACE bnd_getTrac
+  MODULE PROCEDURE bnd_get_trac_
+  MODULE PROCEDURE bnd_get_trac_exc_
+END INTERFACE bnd_getTrac
+
+!> \brief Interface to setter for tractions
+INTERFACE bnd_setTrac
+  MODULE PROCEDURE bnd_set_trac_
+  MODULE PROCEDURE bnd_set_trac_exc_
+END INTERFACE bnd_setTrac
+
+
+
+CONTAINS
+
+
+
+! ************************************************************************
+! ACCESS PROGRAMS
+! ************************************************************************
+
+
+! ------------------------------------------------------------------------
+!>  \brief Constructor for traction element data (non-exception checking)
+!!
+!!  \param  nel  Number of elements
+!!
+!!  This routine allocates memory and initializes the state variable that
+!!  contains the set of elements for the problem.
+! ------------------------------------------------------------------------
+SUBROUTINE bnd_init_ (nel)
+  INTEGER, INTENT(IN) :: nel
+  INTEGER :: i, j     !> loop variables
+
+  !> ensure state variable is clear
+  CALL bnd_clean()
+
+  !> allocate memory for element data state variable
+  ALLOCATE(boundElements(nel))
+
+  !> initialize state variable
+  DO i = 1,nel
+
+    !> element number
+    boundElements(i)%num   = i
+
+    !> initialize connectivity and loads
+    DO j = 1,NNODELB
+      boundElements(i)%ico  %nd(j)            = 0
+      boundElements(i)%trac %loads(j) %sig_nt = 0.d0
+      boundElements(i)%trac %loads(j) %sig_nn = 0.d0
+    END DO
+
+  END DO
+
+END SUBROUTINE bnd_init_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Constructor for traction element data (exception checking)
+!!
+!!  \param  nel   Number of elements
+!!  \param  exc   Error code
+!!
+!!  \exception  ALLOC Failed to allocate memory for element data
+!!  \exception  SZE   Specified number of elements is invalid
+!!
+!!  This routine allocates memory and initializes the state variable that
+!!  contains the set of elements for the problem.
+! ------------------------------------------------------------------------
+SUBROUTINE bnd_init_exc_ (nel, exc)
+  INTEGER, INTENT(IN) :: nel
+  INTEGER, INTENT(OUT) :: exc
+  INTEGER :: e        !> status code for memory allocation
+  INTEGER :: i, j     !> loop variables
+
+  !> ensure state variable is clear
+  CALL bnd_clean()
+
+  !> make sure number of elements is valid
+  IF (nel.LT.1 .OR. nel.GT.MAX_BOUNDELS) THEN
+    exc=SZE
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> allocate memory for element data state variable
+  ALLOCATE(boundElements(nel), STAT=e)
+
+  !> make sure memory allocation was successful
+  IF (e.NE.0) THEN
+    exc=ALLOC
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> initialize state variable
+  DO i = 1,nel
+
+    !> element number
+    boundElements(i)%num   = i
+
+    !> initialize connectivity and loads
+    DO j = 1,NNODELB
+      boundElements(i)%ico  %nd(j)            = 0
+      boundElements(i)%trac %loads(j) %sig_nt = 0.d0
+      boundElements(i)%trac %loads(j) %sig_nn = 0.d0
+    END DO
+
+  END DO
+
+END SUBROUTINE bnd_init_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Destructor for traction element data
+!!
+!!  This routine clears the memory allocated to the state variable that
+!!  contains the set of element data for the problem.
+! ------------------------------------------------------------------------
+SUBROUTINE bnd_clean_ ()
+
+  IF (ALLOCATED(boundElements)) DEALLOCATE(boundElements)
+
+END SUBROUTINE bnd_clean_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for number of traction elements
+!!
+!!  \return nel  Number of elements
+!!
+!!  This routine determines the number of data entries that have been
+!!  allocated for element data. It does not check that the element data
+!!  has been populated (i.e. changed from initial zero values).
+! ------------------------------------------------------------------------
+FUNCTION bnd_num_bound_elem_ () RESULT(nel)
+  INTEGER :: nel
+
+  !> if data is initialized, return number of elements
+  IF (ALLOCATED(boundElements)) THEN
+    nel = SIZE(boundElements)
+  ELSE
+    nel = 0   !> if not initialized, there are no traction elements
+  END IF
+
+END FUNCTION bnd_num_bound_elem_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for connectivity (non-exception checking)
+!!
+!!  \param  i     Element number
+!!  \param  j     Local node index
+!!
+!!  \return n     Global node number
+!!
+!!  This routine determines the global node number associated with local
+!!  node j for element i.
+! ------------------------------------------------------------------------
+FUNCTION bnd_get_connect_ (i,j) RESULT(n)
+  INTEGER, INTENT(IN) :: i,j
+  INTEGER :: n
+
+  n = boundElements(i)%ico%nd(j)
+
+END FUNCTION bnd_get_connect_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for connectivity (exception checking)
+!!
+!!  \param  i     Node number
+!!  \param  j     Local node index
+!!  \param  exc   Error code
+!!
+!!  \return n     Global node number
+!!
+!!  \exception  POSIT The element number is not in [1..bnd_numBoundElem()]
+!!                    or the local node index is not in [1..NNODELB]
+!!
+!!  This routine determines the global node number associated with local
+!!  node j for element i.
+! ------------------------------------------------------------------------
+FUNCTION bnd_get_connect_exc_ (i,j, exc) RESULT(n)
+  INTEGER, INTENT(IN) :: i,j
+  INTEGER, INTENT(OUT) :: exc
+  INTEGER :: n
+
+  !> check that the element number is within the bounds of the elements
+  !! list and that the local node index is within [1..NNODELB]
+  IF ( i.LT.1 .OR. i.GT.bnd_numBoundElem() &
+        .OR. j.LT.1 .OR. j.GT.NNODELB ) THEN
+    exc=POSIT
+    CALL log_printLogMsg(exc,sdr)
+    n = 0
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  n = bnd_getConnect(i,j)
+
+END FUNCTION bnd_get_connect_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Setter for connectivity (non-exception checking)
+!!
+!!  \param  i     Node number
+!!  \param  j     Local node index
+!!  \param  n     Global node number
+!!
+!!  This routine sets the global node number associated with local
+!!  node j for element i.
+! ------------------------------------------------------------------------
+SUBROUTINE bnd_set_connect_ (i,j, n)
+  INTEGER, INTENT(IN) :: i,j
+  INTEGER, INTENT(IN) :: n
+
+  boundElements(i)%ico%nd(j) = n
+
+END SUBROUTINE bnd_set_connect_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Setter for connectivity (exception checking)
+!!
+!!  \param  i     Node number
+!!  \param  j     Local node index
+!!  \param  n     Global node number
+!!  \param  exc   Error code
+!!
+!!  \exception  EXCEED  The value of the node number exceeds the number of
+!!                      nodes
+!!  \exception  POSIT   The element number is not in
+!!                      [1..bnd_numBoundElem()] or the local node index is
+!!                      not in [1..NNODELB]
+!!
+!!  This routine sets the global node number associated with local
+!!  node j for element i.
+! ------------------------------------------------------------------------
+SUBROUTINE bnd_set_connect_exc_ (i,j, n, exc)
+  INTEGER, INTENT(IN) :: i,j
+  INTEGER, INTENT(IN) :: n
+  INTEGER, INTENT(OUT) :: exc
+
+  !> check that the element number is within the bounds of the elements
+  !! list and that the local node index is within [1..NNODELB]
+  IF ( i.LT.1 .OR. i.GT.bnd_numBoundElem() &
+        .OR. j.LT.1 .OR. j.GT.NNODELB ) THEN
+
+    exc=POSIT
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+
+  ELSE IF (n.LT.1 .OR. n.GT.fld_numNode()) THEN
+
+    exc=EXCEED
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  CALL bnd_setConnect(i,j, n)
+
+END SUBROUTINE bnd_set_connect_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Compute the length of an element (non-exception checking)
+!!
+!!  \param  i       Element number
+!!
+!!  \return length  Volume of element
+!!
+!!  This routine computes the length of a traction element. This is
+!!  computed as:
+!!
+!!    L   =  ( (y2-y1)**2 + (x2-x1)**2 )**0.5
+!!
+!!  where xj and yj are the (x,y) coordinates of node j of element i.
+! ------------------------------------------------------------------------
+FUNCTION bnd_len_bound_elem_ (i) RESULT(length)
+  INTEGER, INTENT(IN) :: i
+  DOUBLE PRECISION :: length
+  DOUBLE PRECISION, DIMENSION(NNODELB,NDIM) :: x  !> element coords
+  INTEGER :: j,k    !> loop variable
+
+  !> get element coords
+  DO k = 1,NDIM
+    DO j = 1,NNODELB
+      x(j,k) = fld_getCoord(boundElements(i)%ico%nd(j), k)
+    END DO
+  END DO
+
+  !> compute length
+  length = SQRT( (x(2,2)-x(1,2))**2 + (x(2,1)-x(1,1))**2 )
+
+END FUNCTION bnd_len_bound_elem_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Compute the length of an element (exception checking)
+!!
+!!  \param  i     Element number
+!!  \param  exc   Error code
+!!
+!!  \return length   Length of element
+!!
+!!  \exception  POSIT   The specified element number is not within the
+!!                      range [1..bnd_numBoundElem()]
+!!
+!!  This routine computes the length of a traction element. This is
+!!  computed as:
+!!
+!!    L   =  ( (y2-y1)**2 + (x2-x1)**2 )**0.5
+!!
+!!  where xj and yj are the (x,y) coordinates of node j of element i.
+! ------------------------------------------------------------------------
+FUNCTION bnd_len_bound_elem_exc_ (i, exc) RESULT(length)
+  INTEGER, INTENT(IN) :: i
+  INTEGER, INTENT(OUT) :: exc
+  DOUBLE PRECISION :: length
+
+  !> make sure the element number is valid
+  IF (i.LT.1 .OR. i.GT.bnd_numBoundElem()) THEN
+    exc=POSIT
+    CALL log_printLogMsg(exc,sdr)
+    length = 0.d0
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  length = bnd_lenBoundElem(i)
+
+END FUNCTION bnd_len_bound_elem_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for traction (non-exception checking)
+!!
+!!  \param  i     Element number
+!!  \param  j     Local node index
+!!
+!!  \return load  Surface load
+!!
+!!  This routine determines the surface load associated with local
+!!  node j for element i.
+! ------------------------------------------------------------------------
+FUNCTION bnd_get_trac_ (i,j) RESULT(load)
+  INTEGER, INTENT(IN) :: i,j
+  TYPE(surfLoadT) :: load
+
+  load%sig_nt = boundElements(i)%trac%loads(j)%sig_nt
+  load%sig_nn = boundElements(i)%trac%loads(j)%sig_nn
+
+END FUNCTION bnd_get_trac_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for traction (exception checking)
+!!
+!!  \param  i     Node number
+!!  \param  j     Local node index
+!!  \param  exc   Error code
+!!
+!!  \return load  Surface load
+!!
+!!  \exception  POSIT The element number is not in [1..bnd_numBoundElem()]
+!!                    or the local node index is not in [1..NNODELB]
+!!
+!!  This routine determines the surface load associated with local
+!!  node j for element i.
+! ------------------------------------------------------------------------
+FUNCTION bnd_get_trac_exc_ (i,j, exc) RESULT(load)
+  INTEGER, INTENT(IN) :: i,j
+  INTEGER, INTENT(OUT) :: exc
+  TYPE(surfLoadT) :: load
+
+  !> check that the element number is within the bounds of the elements
+  !! list and that the local node index is within [1..NNODELB]
+  IF ( i.LT.1 .OR. i.GT.bnd_numBoundElem() &
+        .OR. j.LT.1 .OR. j.GT.NNODELB ) THEN
+    exc=POSIT
+    CALL log_printLogMsg(exc,sdr)
+    load%sig_nt = 0.d0
+    load%sig_nn = 0.d0
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  load = bnd_getTrac(i,j)
+
+END FUNCTION bnd_get_trac_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Setter for traction (non-exception checking)
+!!
+!!  \param  i     Node number
+!!  \param  j     Local node index
+!!  \param  load  Surface load
+!!
+!!  This routine sets the surface load associated with local
+!!  node j for element i.
+! ------------------------------------------------------------------------
+SUBROUTINE bnd_set_trac_ (i,j, load)
+  INTEGER, INTENT(IN) :: i,j
+  TYPE(surfLoadT), INTENT(IN) :: load
+
+  boundElements(i)%trac%loads(j)%sig_nt = load%sig_nt
+  boundElements(i)%trac%loads(j)%sig_nn = load%sig_nn
+
+END SUBROUTINE bnd_set_trac_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Setter for traction (exception checking)
+!!
+!!  \param  i     Node number
+!!  \param  j     Local node index
+!!  \param  load  Surface load
+!!  \param  exc   Error code
+!!
+!!  \exception  EXCEED  A value in the set of surface load stress is not
+!!                      within the specified limits
+!!  \exception  POSIT   The element number is not in
+!!                      [1..bnd_numBoundElem()] or the local node index is
+!!                      not in [1..NNODELB]
+!!
+!!  This routine sets the surface load associated with local
+!!  node j for element i.
+! ------------------------------------------------------------------------
+SUBROUTINE bnd_set_trac_exc_ (i,j, load, exc)
+  INTEGER, INTENT(IN) :: i,j
+  TYPE(surfLoadT), INTENT(IN) :: load
+  INTEGER, INTENT(OUT) :: exc
+
+  !> check that the element number is within the bounds of the elements
+  !! list and that the local node index is within [1..NNODELB]
+  IF ( i.LT.1 .OR. i.GT.bnd_numBoundElem() &
+        .OR. j.LT.1 .OR. j.GT.NNODELB ) THEN
+
+    exc=POSIT
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+
+  ELSE IF ( load%sig_nt.LT.SIG_MIN .OR. load%sig_nt.GT.SIG_MAX &
+            .OR. load%sig_nn.LT.SIG_MIN .OR. load%sig_nn.GT.SIG_MAX ) THEN
+
+    exc=EXCEED
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  CALL bnd_setTrac(i,j, load)
+
+END SUBROUTINE bnd_set_trac_exc_
+
+
+
+! ************************************************************************
+! LOCAL FUNCTIONS
+! ************************************************************************
+
+! none
+
+END MODULE boundary_data

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/text-base/boundary_data_test.f90.svn-base

+! ------------------------------------------------------------------------
+!>  \brief Module for testing Boundary Data module
+! ------------------------------------------------------------------------
+MODULE boundary_data_test
+USE fruit                   !> Unit testing framework
+USE system_constants        !> Global constants
+USE log_message_control     !> Printing log/error messages
+USE log_messages            !> Log/error codes
+USE field_data              !> Field Data module
+USE boundary_data           !> Boundary Data module
+
+IMPLICIT NONE
+
+
+CONTAINS
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for OK exception message on allocation of elements
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param nel        Number of elements
+!!  \param expMsg     Expected log message
+!!  \param actMsg     Actual log message
+!!
+!!  This test initializes the elements state variable in the boundary_data
+!!  module and makes sure that the exception message is OK
+!!  (i.e. allocation did not fail).
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_elem_allocation_MSG
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_elem_allocation_MSG'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: nel=100
+  INTEGER, PARAMETER :: expMsg=OK
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log message file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize module and check the exception
+  CALL bnd_init(nel, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate module
+  CALL log_closeLogFile()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_elem_allocation_MSG
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for SZE exception message on allocation of traction elements
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param expMsg     Expected log message
+!!  \param actMsg     Actual log message
+!!
+!!  This test attempts to initialize the elements state variable in the
+!!  boundary data module with invalid size parameters and verifies that
+!!  the correct exception is returned.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_elem_allocation_SZE
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_elem_allocation_SZE'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: expMsg=SZE
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> try to initialize with nel=0
+  CALL bnd_init(0, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> ensure module is reset
+  CALL bnd_clean()
+
+  !> try to initialize with nel=MAX_BOUNDELS+1
+  CALL bnd_init(MAX_BOUNDELS+1, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate module
+  CALL log_closeLogFile()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_elem_allocation_SZE
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for number of elements when element data is not initialized
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param expected     Expected number of elements
+!!  \param actual       Actual number of elements
+!!
+!!  This test makes sure that the number of elements is returned as 0
+!!  when the element data is not initialized
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_num_elem_not_allocated
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_num_elem_not_allocated'
+  INTEGER, PARAMETER :: expected = 0
+  INTEGER :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> check number of materials
+  actual = bnd_numBoundElem()
+  CALL assertEquals(expected, actual)
+
+END SUBROUTINE test_bnd_num_elem_not_allocated
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for number of elements when element data is initialized
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param expected     Expected number of elements
+!!  \param actual       Actual number of elements
+!!
+!!  This test makes sure that the correct number of elements is returned.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_num_elem_allocated
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_num_elem_allocated'
+  INTEGER, PARAMETER :: expected = 100
+  INTEGER :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the module
+  CALL bnd_init(expected)
+
+  !> check number of nodes
+  actual = bnd_numBoundElem()
+  CALL assertEquals(expected, actual)
+
+  !> deallocate the module
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_num_elem_allocated
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for POSIT exception from bnd_getConnect
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param nel        Number of elements
+!!  \param iel        Element number
+!!  \param jnod       Local node index
+!!  \param expMsg     Expected log message
+!!  \param actMsg     Actual log message
+!!  \param n          Dummy variable for get function return
+!!
+!!  This test checks that a POSIT exception is returned when the requested
+!!  location is not inside the element list or the requested local node
+!!  index is not within [1..NNODELB].
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_get_connect_POSIT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_get_connect_POSIT'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: nel=100
+  INTEGER, PARAMETER :: iel = 50
+  INTEGER, PARAMETER :: jnod = 2
+  INTEGER, PARAMETER :: expMsg = POSIT
+  INTEGER :: actMsg
+  INTEGER :: n
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize nodes
+  CALL bnd_init(nel)
+
+  !> try to get beyond last element
+  n = bnd_getConnect(nel+1,jnod, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get before first element
+  n = bnd_getConnect(0,jnod, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get beyond last local node index
+  n = bnd_getConnect(iel,NNODELB+1, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get before first local node index
+  n = bnd_getConnect(iel,0, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize the log file and deallocate the elements list
+  CALL log_closeLogFile()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_get_connect_POSIT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for POSIT exception from bnd_setConnect
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param nnod       Number of nodes
+!!  \param nel        Number of elements
+!!  \param iel        Element number
+!!  \param jnod       Local node index
+!!  \param expMsg     Expected log message
+!!  \param actMsg     Actual log message
+!!  \param n          Dummy variable for set routine input
+!!
+!!  This test checks that a POSIT exception is returned when the requested
+!!  location is not inside the element list or the requested local node
+!!  index is not within [1..NNODELB].
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_set_connect_POSIT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_set_connect_POSIT'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: nnod=100
+  INTEGER, PARAMETER :: nel=100
+  INTEGER, PARAMETER :: iel = 50
+  INTEGER, PARAMETER :: jnod = 2
+  INTEGER, PARAMETER :: expMsg = POSIT
+  INTEGER :: actMsg
+  INTEGER, PARAMETER :: n = 10
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize nodes and elements
+  CALL bnd_init(nel)
+
+  !> try to set beyond last element
+  CALL bnd_setConnect(nel+1,jnod, n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set before first element
+  CALL bnd_setConnect(0,jnod, n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set beyond last local node index
+  CALL bnd_setConnect(iel,NNODELB+1, n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set before first local node index
+  CALL bnd_setConnect(iel,0, n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize the log file and deallocate elements
+  CALL log_closeLogFile()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_set_connect_POSIT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for EXCEED exception from bnd_setConnect
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param nnod       Number of nodes
+!!  \param nel        Number of elements
+!!  \param iel        Element number
+!!  \param jnod       Local node index
+!!  \param expMsg     Expected log message
+!!  \param actMsg     Actual log message
+!!
+!!  This test checks that an EXCEED exception is returned when the
+!!  specified input is not within the range defined in the System
+!!  Constants module.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_set_connect_EXCEED
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_set_connect_EXCEED'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: nnod=100
+  INTEGER, PARAMETER :: nel=100
+  INTEGER, PARAMETER :: iel = 50
+  INTEGER, PARAMETER :: jnod = 2
+  INTEGER, PARAMETER :: expMsg = EXCEED
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize nodes and elements
+  CALL fld_initNode(nnod)
+  CALL bnd_init(nel)
+
+  !> try to set below minimum
+  CALL bnd_setConnect(iel,jnod, 0, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set above maximum
+  CALL bnd_setConnect(iel,jnod, nnod+1, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize the log file and deallocate the nodes and elements
+  CALL log_closeLogFile()
+  CALL fld_cleanNode()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_set_connect_EXCEED
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct value getting and setting in bnd_getConnect
+!!        and bnd_setConnect
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param nnod       Number of nodes
+!!  \param nel        Number of elements
+!!  \param iel        Element number
+!!  \param jnod       Local node index
+!!  \param expVal     Expected data value
+!!  \param actVal     Actual data value
+!!
+!!  This test checks that the correct value is set using the
+!!  bnd_setConnect access program and returned from the
+!!  bnd_getConnect access program.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_get_set_connect_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_get_set_connect_VAL'
+  INTEGER, PARAMETER :: nnod=100
+  INTEGER, PARAMETER :: nel=100
+  INTEGER, PARAMETER :: iel = 50
+  INTEGER, PARAMETER :: jnod = 2
+  INTEGER, PARAMETER :: expVal = 20
+  INTEGER :: actVal
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize nodes and elements
+  CALL fld_initNode(nnod)
+  CALL bnd_init(nel)
+
+  !> set connectivity value
+  CALL bnd_setConnect(iel,jnod, expVal)
+  actVal = bnd_getConnect(iel,jnod)
+  CALL assertEquals(expVal,actVal)
+
+  !> finalize the log file and deallocate the nodes and elements
+  CALL log_closeLogFile()
+  CALL fld_cleanNode()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_get_set_connect_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for POSIT exception message from bnd_lenBoundElem
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param nel        Number of elements
+!!  \param expMsg     Expected log message
+!!  \param actMsg     Actual log message
+!!  \param length     Dummy variable for get function return
+!!
+!!  This test makes sure that the correct exception is returned when the
+!!  specified element index is not in the range of the allocated element
+!!  data.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_len_bound_elem_POSIT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_len_bound_elem_POSIT'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: nel=100
+  INTEGER, PARAMETER :: expMsg=POSIT
+  INTEGER :: actMsg
+  DOUBLE PRECISION :: length
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize elements
+  CALL bnd_init(nel)
+
+  !> try to get beyond last element
+  length = bnd_lenBoundElem(nel+1, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get before first element
+  length = bnd_lenBoundElem(0, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize the log file and deallocate the elements list
+  CALL log_closeLogFile()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_len_bound_elem_POSIT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for element length when node coordinates not initialized
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param nnod       Number of nodes
+!!  \param nel        Number of elements
+!!  \param iel        Element number
+!!  \param expected   Expected length of element (zero)
+!!  \param actual     Actual length of element
+!!
+!!  This test makes sure that the length of the element is returned as 0
+!!  when the coordinates of the connected nodes have not been set
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_len_bound_elem_ZERO
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_len_bound_elem_ZERO'
+  INTEGER, PARAMETER :: nnod=100
+  INTEGER, PARAMETER :: nel=100
+  DOUBLE PRECISION, PARAMETER :: expected = 0.d0
+  DOUBLE PRECISION :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize nodes and elements
+  CALL fld_initNode(nnod)
+  CALL bnd_init(nel)
+
+  !> set connectivity of first element to {1,2}
+  CALL bnd_setConnect(1,1,1)
+  CALL bnd_setConnect(1,2,2)
+
+  !> check length of element
+  actual = bnd_lenBoundElem(1)
+  CALL assertEquals(expected,actual)
+
+  !> deallocate the elements and nodes
+  CALL fld_cleanNode()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_len_bound_elem_ZERO
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for element length when node coordinates initialized
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param nnod       Number of nodes
+!!  \param nel        Number of elements
+!!  \param iel        Element number
+!!  \param expected   Expected length of element
+!!  \param actual     Actual length of element
+!!
+!!  This test makes sure that the length of the element is computed
+!!  correctly for the line given by (1,1), (3,2)
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_len_bound_elem_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_len_bound_elem_VAL'
+  INTEGER, PARAMETER :: nnod=100
+  INTEGER, PARAMETER :: nel=100
+  DOUBLE PRECISION, PARAMETER :: toler = 1.d-14
+  DOUBLE PRECISION, PARAMETER :: expected = 2.23606797749979d0
+  DOUBLE PRECISION :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize nodes and elements
+  CALL fld_initNode(nnod)
+  CALL bnd_init(nel)
+
+  !> set connectivity of first element to {1,2}
+  CALL bnd_setConnect(1,1,1)
+  CALL bnd_setConnect(1,2,2)
+
+  !> set coordinates of nodes
+  CALL fld_setCoord(1,1, 1.d0)
+  CALL fld_setCoord(1,2, 1.d0)
+  CALL fld_setCoord(2,1, 3.d0)
+  CALL fld_setCoord(2,2, 2.d0)
+
+  !> check volume of element
+  actual = bnd_lenBoundElem(1)
+  CALL assertEquals(expected,actual, toler)
+
+  !> deallocate the elements and nodes
+  CALL fld_cleanNode()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_len_bound_elem_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for POSIT exception from bnd_getTrac
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param nel        Number of elements
+!!  \param iel        Element number
+!!  \param jnod       Local node index
+!!  \param expMsg     Expected log message
+!!  \param actMsg     Actual log message
+!!  \param load       Dummy variable for get function return
+!!
+!!  This test checks that a POSIT exception is returned when the requested
+!!  location is not inside the element list or the requested local node
+!!  index is not within [1..NNODELB].
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_get_trac_POSIT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_get_trac_POSIT'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: nel=100
+  INTEGER, PARAMETER :: iel = 50
+  INTEGER, PARAMETER :: jnod = 2
+  INTEGER, PARAMETER :: expMsg = POSIT
+  INTEGER :: actMsg
+  TYPE(surfLoadT) :: load
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize nodes
+  CALL bnd_init(nel)
+
+  !> try to get beyond last element
+  load = bnd_getTrac(nel+1,jnod, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get before first element
+  load = bnd_getTrac(0,jnod, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get beyond last local node index
+  load = bnd_getTrac(iel,NNODELB+1, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get before first local node index
+  load = bnd_getTrac(iel,0, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize the log file and deallocate the elements list
+  CALL log_closeLogFile()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_get_trac_POSIT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for POSIT exception from bnd_setTrac
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param nnod       Number of nodes
+!!  \param nel        Number of elements
+!!  \param iel        Element number
+!!  \param jnod       Local node index
+!!  \param expMsg     Expected log message
+!!  \param actMsg     Actual log message
+!!  \param load       Dummy variable for set routine input
+!!
+!!  This test checks that a POSIT exception is returned when the requested
+!!  location is not inside the element list or the requested local node
+!!  index is not within [1..NNODELB].
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_set_trac_POSIT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_set_trac_POSIT'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: nnod=100
+  INTEGER, PARAMETER :: nel=100
+  INTEGER, PARAMETER :: iel = 50
+  INTEGER, PARAMETER :: jnod = 2
+  INTEGER, PARAMETER :: expMsg = POSIT
+  INTEGER :: actMsg
+  TYPE(surfLoadT) :: load
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize elements
+  CALL bnd_init(nel)
+
+  !> try to set beyond last element
+  CALL bnd_setTrac(nel+1,jnod, load, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set before first element
+  CALL bnd_setTrac(0,jnod, load, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set beyond last local node index
+  CALL bnd_setTrac(iel,NNODELB+1, load, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set before first local node index
+  CALL bnd_setTrac(iel,0, load, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize the log file and deallocate the nodes and elements
+  CALL log_closeLogFile()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_set_trac_POSIT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for EXCEED exception from bnd_setTrac
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param nel        Number of elements
+!!  \param iel        Element number
+!!  \param jnod       Local node index
+!!  \param sig_nt     Test shear stress
+!!  \param sig_nn     Test normal stress
+!!  \param load       Surface load
+!!  \param expMsg     Expected log message
+!!  \param actMsg     Actual log message
+!!
+!!  This test checks that an EXCEED exception is returned when the
+!!  specified input is not within the range defined in the System
+!!  Constants module.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_set_trac_EXCEED
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_set_trac_EXCEED'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: nel=100
+  INTEGER, PARAMETER :: iel = 50
+  INTEGER, PARAMETER :: jnod = 2
+  DOUBLE PRECISION :: sig_nt, sig_nn
+  TYPE(surfLoadT) :: load
+  INTEGER, PARAMETER :: expMsg = EXCEED
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize elements
+  CALL bnd_init(nel)
+
+  !> set load such that sig_nt < SIG_MIN
+  sig_nt = SIG_MIN-0.1d0
+  sig_nn = 0.d0
+  load%sig_nt = sig_nt
+  load%sig_nn = sig_nn
+
+  !> try to set traction
+  CALL bnd_setTrac(iel,jnod, load, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> set load such that sig_nt > SIG_MAX
+  sig_nt = SIG_MAX+0.1d0
+  load%sig_nt = sig_nt
+
+  !> try to set traction
+  CALL bnd_setTrac(iel,jnod, load, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> set load such that sig_nn < SIG_MIN
+  sig_nt = 0.d0
+  sig_nn = SIG_MIN-0.1d0
+  load%sig_nt = sig_nt
+  load%sig_nn = sig_nn
+
+  !> try to set traction
+  CALL bnd_setTrac(iel,jnod, load, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> set load such that sig_nn > SIG_MAX
+  sig_nn = SIG_MAX+0.1d0
+  load%sig_nn = sig_nn
+
+  !> try to set traction
+  CALL bnd_setTrac(iel,jnod, load, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize the log file and deallocate elements
+  CALL log_closeLogFile()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_set_trac_EXCEED
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct value getting and setting in bnd_getTrac
+!!        and bnd_setTrac
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param nel        Number of elements
+!!  \param iel        Element number
+!!  \param jnod       Local node index
+!!  \param load       Surface load
+!!  \param expSigNT   Expected shear stress
+!!  \param expSigNN   Expected normal stress
+!!  \param actSigNT   Actual shear stress
+!!  \param actSigNN   Actual normal stress
+!!
+!!  This test checks that the correct value is set using the
+!!  bnd_setTrac access program and returned from the
+!!  bnd_getTrac access program.
+! ------------------------------------------------------------------------
+SUBROUTINE test_bnd_get_set_trac_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_bnd_get_set_trac_VAL'
+  INTEGER, PARAMETER :: nnod=100
+  INTEGER, PARAMETER :: nel=100
+  INTEGER, PARAMETER :: iel = 50
+  INTEGER, PARAMETER :: jnod = 2
+  DOUBLE PRECISION, PARAMETER :: expSigNT=1.d3, expSigNN=-2.d3
+  TYPE(surfLoadT) :: load
+  DOUBLE PRECISION :: actSigNT, actSigNN
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize elements
+  CALL bnd_init(nel)
+
+  !> set traction value
+  load%sig_nt = expSigNT
+  load%sig_nn = expSigNN
+  CALL bnd_setTrac(iel,jnod, load)
+
+  !> reset load
+  load%sig_nt = 0.d0
+  load%sig_nn = 0.d0
+
+  !> get traction value
+  load = bnd_getTrac(iel,jnod)
+  actSigNT = load%sig_nt
+  actSigNN = load%sig_nn
+  CALL assertEquals(expSigNT,actSigNT)
+  CALL assertEquals(expSigNN,actSigNN)
+
+  !> finalize the log file and deallocate the nodes and elements
+  CALL log_closeLogFile()
+  CALL bnd_clean()
+
+END SUBROUTINE test_bnd_get_set_trac_VAL
+
+END MODULE boundary_data_test

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/text-base/constitutive.f90.svn-base

+! ------------------------------------------------------------------------
+!>  \brief Module for Constitutive Matrix
+! ------------------------------------------------------------------------
+MODULE constitutive
+USE system_constants    !> Global constants (for size of matrix)
+USE log_message_control !> Printing log/error messages
+USE log_messages        !> Log/error codes
+USE dense_matrix_def    !> Dense Matrix data type
+
+IMPLICIT NONE
+PRIVATE
+
+
+
+! ************************************************************************
+! EXPORTS
+! ************************************************************************
+
+
+!> Exported interfaces
+PUBLIC :: dmatrix
+
+
+
+! ************************************************************************
+! LOCAL CONSTANTS
+! ************************************************************************
+
+INTEGER, PARAMETER :: sdr=CNSMAT
+
+
+! ************************************************************************
+! DATA TYPES
+! ************************************************************************
+
+! none
+
+
+
+! ************************************************************************
+! INTERFACES
+! ************************************************************************
+
+
+!> \brief Interface to constitutive matrix access program
+INTERFACE dmatrix
+  MODULE PROCEDURE dmatrix_linear_elastic_plane_strain_
+  MODULE PROCEDURE dmatrix_linear_elastic_plane_strain_exc_
+END INTERFACE dmatrix
+
+
+
+CONTAINS
+
+
+
+! ************************************************************************
+! ACCESS PROGRAMS
+! ************************************************************************
+
+
+! ------------------------------------------------------------------------
+!>  \brief Constitutive matrix (linear elastic, plane strain)
+!!
+!!  \param  emod  Elastic modulus
+!!  \param  nu    Poisson's ratio
+!!  \param  Dmat  Constitutive matrix
+!!
+!!  This routine builds the constitutive matrix for a linear elastic
+!!  material under plane strain conditions. The matrix has the following
+!!  form:
+!!
+!!    D = ( emod / ((1+nu)*(1-2*nu)) ) * [ 1-nu   nu        0       ]
+!!                                       [  nu   1-nu       0       ]
+!!                                       [   0     0   0.5*(1-2*v)  ]
+! ------------------------------------------------------------------------
+SUBROUTINE dmatrix_linear_elastic_plane_strain_ (emod,nu, Dmat)
+  DOUBLE PRECISION, INTENT(IN) :: emod,nu
+  TYPE(matrixT), INTENT(INOUT) :: Dmat
+  DOUBLE PRECISION :: one_minus_nu, one_minus_two_nu, coef
+
+  one_minus_nu = 1.d0 - nu
+  one_minus_two_nu = 1.d0 - 2.d0*nu
+  coef = emod / ((1.d0+nu)*one_minus_two_nu)
+
+  CALL dm_init(Dmat, NTNS,NTNS)
+
+  CALL dm_set(Dmat, 1,1, coef*one_minus_nu)
+  CALL dm_set(Dmat, 2,1, coef*nu)
+  CALL dm_set(Dmat, 1,2, dm_get(Dmat, 2,1))
+  CALL dm_set(Dmat, 2,2, dm_get(Dmat, 1,1))
+  CALL dm_set(Dmat, 3,3, 0.5d0*coef*one_minus_two_nu)
+
+END SUBROUTINE dmatrix_linear_elastic_plane_strain_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Constitutive matrix (linear elastic, plane strain)
+!!
+!!  \param  emod  Elastic modulus
+!!  \param  nu    Poisson's ratio
+!!  \param  Dmat  Constitutive matrix
+!!  \param  exc   Error code
+!!
+!!  \exception  EXCEED  The specified values of elastic modulus or
+!!                      Poisson's ratio are not within specified limits
+!!
+!!  This routine builds the constitutive matrix for a linear elastic
+!!  material under plane strain conditions. The matrix has the following
+!!  form:
+!!
+!!    D = ( emod / ((1+nu)*(1-2*nu)) ) * [ 1-nu   nu        0       ]
+!!                                       [  nu   1-nu       0       ]
+!!                                       [   0     0   0.5*(1-2*v)  ]
+! ------------------------------------------------------------------------
+SUBROUTINE dmatrix_linear_elastic_plane_strain_exc_ (emod,nu, Dmat, exc)
+  DOUBLE PRECISION, INTENT(IN) :: emod,nu
+  TYPE(matrixT), INTENT(INOUT) :: Dmat
+  INTEGER, INTENT(OUT) :: exc
+  DOUBLE PRECISION :: one_minus_nu, one_minus_two_nu, coef
+
+  !> check that specified parameters are within limits
+  IF (  emod.LT.E_MIN .OR. emod.GT.E_MAX  &
+        .OR. nu.LT.NU_MIN .OR. nu.GT.NU_MAX ) THEN
+    exc=EXCEED
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  CALL dmatrix(emod,nu, Dmat)
+
+END SUBROUTINE dmatrix_linear_elastic_plane_strain_exc_
+
+
+
+! ************************************************************************
+! LOCAL FUNCTIONS
+! ************************************************************************
+
+! none
+
+END MODULE constitutive

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/text-base/constitutive_test.f90.svn-base

+! ------------------------------------------------------------------------
+!>  \brief Module for testing Constitutive Matrix module
+! ------------------------------------------------------------------------
+MODULE constitutive_test
+USE fruit                   !> Unit testing framework
+USE system_constants        !> Global system constants
+USE log_message_control     !> Printing log/error messages
+USE log_messages            !> Log/error message and sender codes
+USE dense_matrix_def        !> Dense Matrix ADT
+USE constitutive            !> Constitutive Matrix module
+
+IMPLICIT NONE
+
+
+CONTAINS
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for EXCEED exception from dmatrix
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param emod       Non-exceeding elastic modulus
+!!  \param nu         Non-exceeding Poisson's ratio
+!!  \param emod1      Test elastic modulus 1 ( = E_MIN-0.1 < E_MIN )
+!!  \param nu1        Test Poisson's ratio 1 ( = NU_MIN-0.1 < NU_MIN )
+!!  \param emod1      Test elastic modulus 2 ( = EMAX+1.0 > E_MAX )
+!!  \param nu1        Test Poisson's ratio 2 ( = 0.5 > NU_MAX )
+!!  \param expMsg     Expected log message
+!!  \param actMsg     Actual log message
+!!  \param Dmat       Dummy matrix for dmatrix output
+!!
+!!  This test attempts to obtain the constitutive matrix for values that
+!!  are outside the specified range and ensures that the correct exception
+!!  is raised.
+! ------------------------------------------------------------------------
+SUBROUTINE test_constitutive_EXCEED
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_constitutive_EXCEED'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  DOUBLE PRECISION, PARAMETER :: emod = 70.d3, nu = 0.3d0
+  DOUBLE PRECISION, PARAMETER :: emod1 = E_MIN-0.1d0
+  DOUBLE PRECISION, PARAMETER :: nu1 = NU_MIN-0.1d0
+  DOUBLE PRECISION, PARAMETER :: emod2 = E_MAX+1.d0
+  DOUBLE PRECISION, PARAMETER :: nu2 = 0.5d0
+  INTEGER, PARAMETER :: expMsg=EXCEED
+  INTEGER :: actMsg
+  TYPE(matrixT) :: Dmat
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log message file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> try to call with low elastic modulus
+  CALL dmatrix(emod1,nu, Dmat, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to call with high elastic modulus
+  CALL dmatrix(emod2,nu, Dmat, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to call with low Poisson's ratio
+  CALL dmatrix(emod,nu1, Dmat, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to call with high Poisson's ratio
+  CALL dmatrix(emod,nu2, Dmat, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate D matrix
+  CALL log_closeLogFile()
+  CALL dm_clean(Dmat)
+
+END SUBROUTINE test_constitutive_EXCEED
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for dmatrix when elastic modulus is zero
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param emod       Elastic modulus ( = 0 )
+!!  \param nu         Poisson's ratio
+!!  \param expDat     Expected constitutive matrix (zeros)
+!!  \param Dmat       Matrix for dmatrix output
+!!
+!!  This test obtains the constitutive matrix when elastic modulus is
+!!  zero and Poisson's ratio is any value within the system defined
+!!  limits. It should be all zeros.
+! ------------------------------------------------------------------------
+SUBROUTINE test_constitutive_ZERO
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_constitutive_ZERO'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  DOUBLE PRECISION, PARAMETER :: emod = 0.d0
+  DOUBLE PRECISION, PARAMETER :: nu = 0.3d0
+  DOUBLE PRECISION, DIMENSION(NTNS,NTNS) :: expDat
+  TYPE(matrixT) :: Dmat
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result
+  expDat = RESHAPE( (/ &
+                    0.d0, 0.d0, 0.d0, &
+                    0.d0, 0.d0, 0.d0, &
+                    0.d0, 0.d0, 0.d0  &
+                    /), SHAPE(expDat) )
+
+  !> get constitutive matrix
+  CALL dmatrix(emod,nu, Dmat)
+  CALL assertEquals(expDat, Dmat%dat, NTNS,NTNS)
+
+  !> deallocate constitutive matrix
+  CALL dm_clean(Dmat)
+
+END SUBROUTINE test_constitutive_ZERO
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for dmatrix when parameters are in normal range
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param emod       Elastic modulus
+!!  \param nu         Poisson's ratio
+!!  \param toler      Tolerance on output
+!!  \param expDat     Expected constitutive matrix
+!!  \param Dmat       Matrix for dmatrix output
+!!
+!!  This test obtains the constitutive matrix when elastic modulus
+!!  and Poisson's ratio are given values within the system defined
+!!  limits.
+! ------------------------------------------------------------------------
+SUBROUTINE test_constitutive_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_constitutive_VAL'
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  DOUBLE PRECISION, PARAMETER :: emod = 70.d3
+  DOUBLE PRECISION, PARAMETER :: nu = 0.3d0
+  DOUBLE PRECISION, PARAMETER :: toler = 1.d-9
+  DOUBLE PRECISION, DIMENSION(NTNS,NTNS) :: expDat
+  TYPE(matrixT) :: Dmat
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result
+  expDat = RESHAPE( (/ &
+                    94230.7692307692d0, 40384.6153846154d0, 0.d0, &
+                    40384.6153846154d0, 94230.7692307692d0, 0.d0, &
+                    0.d0, 0.d0, 26923.0769230769d0 &
+                    /), SHAPE(expDat) )
+
+  !> get constitutive matrix
+  CALL dmatrix(emod,nu, Dmat)
+  CALL assertEquals(expDat, Dmat%dat, NTNS,NTNS, toler)
+
+  !> deallocate constitutive matrix
+  CALL dm_clean(Dmat)
+
+END SUBROUTINE test_constitutive_VAL
+
+END MODULE constitutive_test

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/text-base/dense_matrix_def.f90.svn-base

+! ------------------------------------------------------------------------
+!>  \brief Module defining Dense Matrix data type
+! ------------------------------------------------------------------------
+MODULE dense_matrix_def
+USE log_message_control !> Print log/error messages
+USE log_messages        !> Log/error codes
+USE vector_def          !> Vector data type
+
+IMPLICIT NONE
+PRIVATE
+
+
+
+! ************************************************************************
+! EXPORTS
+! ************************************************************************
+
+
+!> Exported data types
+PUBLIC :: matrixT
+
+!> Exported interfaces
+PUBLIC :: dm_init, dm_clean, &
+          dm_numRows, dm_numCols, &
+          dm_get, dm_set, &
+          dm_add, OPERATOR (+), &
+          dm_scalMul, dm_vecMul, dm_matMul, OPERATOR (*), &
+          dm_transpose
+
+
+
+! ************************************************************************
+! LOCAL CONSTANTS
+! ************************************************************************
+
+
+!> Sender code for DeNSe MATrix module
+INTEGER, PARAMETER :: sdr = DNSMAT
+
+
+
+! ************************************************************************
+! DATA TYPES
+! ************************************************************************
+
+
+! ------------------------------------------------------------------------
+!>  \brief Dense Matrix ADT structure
+!!
+!!  \param dat  Array to store matrix data
+! ------------------------------------------------------------------------
+TYPE matrixT
+  DOUBLE PRECISION, ALLOCATABLE :: dat(:,:)
+END TYPE matrixT
+
+
+
+! ************************************************************************
+! INTERFACES
+! ************************************************************************
+
+
+!>  \brief Interface to constructor
+INTERFACE dm_init
+  MODULE PROCEDURE dm_init_
+  MODULE PROCEDURE dm_init_exc_
+END INTERFACE dm_init
+
+!>  \brief Interface to destructor
+INTERFACE dm_clean
+  MODULE PROCEDURE dm_clean_
+END INTERFACE dm_clean
+
+!>  \brief Interface to number of rows
+INTERFACE dm_numRows
+  MODULE PROCEDURE dm_num_rows_
+END INTERFACE dm_numRows
+
+!>  \brief Interface to number of columns
+INTERFACE dm_numCols
+  MODULE PROCEDURE dm_num_cols_
+END INTERFACE dm_numCols
+
+!>  \brief Interface to getter for individual entries in matrix
+INTERFACE dm_get
+  MODULE PROCEDURE dm_get_
+  MODULE PROCEDURE dm_get_exc_
+END INTERFACE dm_get
+
+!>  \brief Interface to setter for individual entries in matrix
+INTERFACE dm_set
+  MODULE PROCEDURE dm_set_
+  MODULE PROCEDURE dm_set_exc_
+END INTERFACE dm_set
+
+!>  \brief Interface for addition of two matrixT
+INTERFACE dm_add
+  MODULE PROCEDURE dm_add_
+  MODULE PROCEDURE dm_add_exc_
+END INTERFACE dm_add
+
+!>  \brief Operator overload for addition
+INTERFACE OPERATOR (+)
+  MODULE PROCEDURE dm_add_
+END INTERFACE
+
+!>  \brief Interface for scalar multiplication
+INTERFACE dm_scalMul
+  MODULE PROCEDURE dm_scal_mul_mat_scal_
+END INTERFACE dm_scalMul
+
+!>  \brief Interface for matrix-vector multiplication
+INTERFACE dm_vecMul
+  MODULE PROCEDURE dm_vec_mul_
+  MODULE PROCEDURE dm_vec_mul_exc_
+END INTERFACE dm_vecMul
+
+!>  \brief Interface for matrix-matrix multiplication
+INTERFACE dm_matMul
+  MODULE PROCEDURE dm_mat_mul_
+  MODULE PROCEDURE dm_mat_mul_exc_
+END INTERFACE dm_matMul
+
+!>  \brief Operator overload for multiplication
+!!  (includes scalar, matrix-vector, and matrix-matrix multiplication)
+INTERFACE OPERATOR (*)
+  MODULE PROCEDURE dm_scal_mul_mat_scal_
+  MODULE PROCEDURE dm_scal_mul_scal_mat_
+  MODULE PROCEDURE dm_vec_mul_
+  MODULE PROCEDURE dm_mat_mul_
+END INTERFACE
+
+!>  \brief Interface for matrix transpose
+INTERFACE dm_transpose
+  MODULE PROCEDURE dm_transpose_
+END INTERFACE dm_transpose
+
+
+
+CONTAINS
+
+
+
+! ************************************************************************
+! ACCESS PROGRAMS
+! ************************************************************************
+
+
+! ------------------------------------------------------------------------
+!>  \brief Constructor for matrixT (non-exception checking)
+!!
+!!  \param self Reference to the dense matrix to be initialized
+!!  \param m    Number of rows in the matrix
+!!  \param n    Number of columns in the matrix
+!!
+!!  This routine initializes the matrixT object referenced by self.
+!!  The data structures contained in self are allocated and initial values
+!!  are set to zero.
+! ------------------------------------------------------------------------
+SUBROUTINE dm_init_ (self, m,n)
+  TYPE(matrixT), INTENT(INOUT) :: self
+  INTEGER, INTENT(IN) :: m,n
+  INTEGER :: i,j    !> loop variables
+
+  !> only reallocate if new dimensions do not match existing dimensions
+  IF (dm_numRows(self).NE.m .OR. dm_numCols(self).NE.n) THEN
+
+    !> reset object if it was previously initialized
+    CALL dm_clean(self)
+
+    !> allocate new data structure
+    ALLOCATE(self%dat(m,n))
+
+  END IF
+
+  !> ensure data is initialized to zero
+  DO j = 1,n
+    DO i = 1,m
+      self%dat(i,j) = 0.d0
+    END DO
+  END DO
+
+END SUBROUTINE dm_init_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Constructor for matrixT (exception checking)
+!!
+!!  \param self Reference to the dense matrix to be initialized
+!!  \param m    Number of rows in the matrix
+!!  \param n    Number of columns in the matrix
+!!  \param exc  Error code
+!!
+!!  This routine initializes the matrixT object referenced by self.
+!!  The data structures contained in self are allocated and initial values
+!!  are set to zero. If the specified dimensions are invalid, a SZE
+!!  exception is returned. If allocation of data structures fails, an
+!!  ALLOC exception is returned.
+! ------------------------------------------------------------------------
+SUBROUTINE dm_init_exc_ (self, m,n, exc)
+  TYPE(matrixT), INTENT(INOUT) :: self
+  INTEGER, INTENT(IN) :: m,n
+  INTEGER, INTENT(OUT) :: exc
+  INTEGER :: e      !> allocation info code
+  INTEGER :: i,j    !> loop variables
+
+  !> if dimensions are invalid, raise SZE exception
+  IF (m.LE.0 .OR. n.LE.0) THEN
+    exc=SZE
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> only reallocate if new dimensions do not match existing dimensions
+  IF (dm_numRows(self).NE.m .OR. dm_numCols(self).NE.n) THEN
+
+    !> reset object if it was previously initialized
+    CALL dm_clean(self)
+
+    !> allocate new data structure for matrix storage
+    ALLOCATE(self%dat(m,n), STAT=e)
+
+    !> if error code is returned, raise ALLOC exception
+    IF (e.NE.0) THEN
+      exc=ALLOC
+      CALL log_printLogMsg(exc,sdr)
+      RETURN
+    ELSE
+      exc=OK
+    END IF
+
+  END IF
+
+  !> ensure data is initialized to zero
+  DO j = 1,n
+    DO i = 1,m
+      self%dat(i,j) = 0.d0
+    END DO
+  END DO
+
+END SUBROUTINE dm_init_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Destructor for matrixT
+!!
+!!  \param self Reference to the dense matrix to be destroyed
+!!
+!!  This routine deallocates existing data structures.
+! ------------------------------------------------------------------------
+SUBROUTINE dm_clean_ (self)
+  TYPE(matrixT), INTENT(INOUT) :: self
+
+  IF (ALLOCATED(self%dat)) DEALLOCATE(self%dat)
+
+END SUBROUTINE dm_clean_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for number of rows
+!!
+!!  \param  self  Reference to the dense matrix object
+!!  \return m     Number of rows in the matrix
+!!
+!!  This routine determines the number of rows allocated to the matrix
+!!  object. If the matrix is not initialized it returns 0.
+! ------------------------------------------------------------------------
+FUNCTION dm_num_rows_ (self) RESULT(m)
+  TYPE(matrixT), INTENT(IN) :: self
+  INTEGER :: m
+
+  !> if data is initialized, return number of rows
+  IF (ALLOCATED(self%dat)) THEN
+    m = SIZE(self%dat,1)
+  ELSE
+    m = 0   !> if not initialized, there are no rows
+  END IF
+
+END FUNCTION dm_num_rows_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for number of columns
+!!
+!!  \param  self  Reference to the dense matrix object
+!!  \return n     Number of columns in the matrix
+!!
+!!  This routine determines the number of columns allocated to the matrix
+!!  object. If the matrix is not initialized it returns 0.
+! ------------------------------------------------------------------------
+FUNCTION dm_num_cols_ (self) RESULT(n)
+  TYPE(matrixT), INTENT(IN) :: self
+  INTEGER :: n
+
+  !> if data is initialized, return number of columns
+  IF (ALLOCATED(self%dat)) THEN
+    n = SIZE(self%dat,2)
+  ELSE
+    n = 0   !> if not initialized, there are no columns
+  END IF
+
+END FUNCTION dm_num_cols_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for matrix entries (non-exception checking)
+!!
+!!  \param  self  Reference to the dense matrix object
+!!  \param  i     Row index
+!!  \param  j     Column index
+!!  \return v     Value at location (i,j) of the matrix
+!!
+!!  This routine determines the value at a particular location in the
+!!  matrix.
+! ------------------------------------------------------------------------
+FUNCTION dm_get_ (self, i,j) RESULT(v)
+  TYPE(matrixT), INTENT(IN) :: self
+  INTEGER, INTENT(IN) :: i,j
+  DOUBLE PRECISION :: v
+
+  v = self%dat(i,j)
+
+END FUNCTION dm_get_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Getter for matrix entries (exception checking)
+!!
+!!  \param  self  Reference to the dense matrix object
+!!  \param  i     Row index
+!!  \param  j     Column index
+!!  \param  exc   Error code
+!!  \return v     Value at location (i,j) of the matrix
+!!
+!!  This routine determines the value at a particular location in the
+!!  matrix. If the location is not inside the matrix, it returns a POSIT
+!!  exception.
+! ------------------------------------------------------------------------
+FUNCTION dm_get_exc_ (self, i,j, exc) RESULT(v)
+  TYPE(matrixT), INTENT(IN) :: self
+  INTEGER, INTENT(IN) :: i,j
+  INTEGER, INTENT(OUT) :: exc
+  DOUBLE PRECISION :: v
+
+  !> make sure desired indices are inside the matrix
+  IF (i.GT.dm_numRows(self) .OR. i.LE.0 &
+      .OR. j.GT.dm_numCols(self) .OR. j.LE.0) THEN
+    exc=POSIT
+    CALL log_printLogMsg(exc,sdr)
+    v = 0.d0
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception getter
+  v = dm_get(self, i,j)
+
+END FUNCTION dm_get_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Setter for matrix entries (non-exception checking)
+!!
+!!  \param self Reference to the dense matrix object
+!!  \param i    Row index
+!!  \param j    Column index
+!!  \param v    Value to be stored at location (i,j) of the matrix
+!!
+!!  This routine sets the value at a particular location in the
+!!  matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE dm_set_ (self, i,j,v)
+  TYPE(matrixT), INTENT(INOUT) :: self
+  INTEGER, INTENT(IN) :: i,j
+  DOUBLE PRECISION, INTENT(IN) :: v
+
+  self%dat(i,j) = v
+
+END SUBROUTINE dm_set_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Setter for matrix entries (exception checking)
+!!
+!!  \param self Reference to the dense matrix object
+!!  \param i    Row index
+!!  \param j    Column index
+!!  \param v    Value to be stored at location (i,j) of the matrix
+!!  \param exc  Error code
+!!
+!!  This routine sets the value at a particular location in the
+!!  matrix. If the location is not inside the matrix, it returns a POSIT
+!!  exception.
+! ------------------------------------------------------------------------
+SUBROUTINE dm_set_exc_ (self, i,j,v, exc)
+  TYPE(matrixT), INTENT(INOUT) :: self
+  INTEGER, INTENT(IN) :: i,j
+  DOUBLE PRECISION, INTENT(IN) :: v
+  INTEGER, INTENT(OUT) :: exc
+
+  !> make sure desired indices are inside the matrix
+  IF (i.GT.dm_numRows(self) .OR. i.LE.0 &
+      .OR. j.GT.dm_numCols(self) .OR. j.LE.0) THEN
+    exc=POSIT
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception setter
+  CALL dm_set(self, i,j,v)
+
+END SUBROUTINE dm_set_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Add two dense matrix objects (non-exception checking)
+!!
+!!  \param  self  Reference to the first (left-hand) dense matrix object
+!!  \param  other Reference to the second (right-hand) dense matrix object
+!!  \return new   Reference to the resulting banded dense object
+!!
+!!  This routine adds two dense matrix objects.
+! ------------------------------------------------------------------------
+FUNCTION dm_add_ (self,other) RESULT(new)
+  TYPE(matrixT), INTENT(IN) :: self,other
+  TYPE(matrixT) :: new
+
+  !> initialize solution matrix
+  CALL dm_init(new, dm_numRows(self),dm_numCols(self))
+
+  !> Add the matrices (using array ops for efficiency)
+  new%dat = self%dat + other%dat
+
+END FUNCTION dm_add_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Add two dense matrix objects (exception checking)
+!!
+!!  \param  self  Reference to the first (left-hand) dense matrix object
+!!  \param  other Reference to the second (right-hand) dense matrix object
+!!  \param  exc   Error code
+!!  \return new   Reference to the resulting dense matrix object
+!!
+!!  This routine adds two dense matrix objects. If the dimensions of the
+!!  two input matrices do not match, it returns a DIMEN exception.
+! ------------------------------------------------------------------------
+FUNCTION dm_add_exc_ (self,other, exc) RESULT(new)
+  TYPE(matrixT), INTENT(IN) :: self,other
+  INTEGER, INTENT(OUT) :: exc
+  TYPE(matrixT) :: new
+
+  !> make sure that the two matrices have the same number of rows
+  IF(dm_numRows(self).NE.dm_numRows(other) &
+          .OR. dm_numCols(self).NE.dm_numCols(other)) THEN
+    exc=DIMEN
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  new = dm_add(self,other)
+
+END FUNCTION dm_add_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Scalar multiplication (self*k version)
+!!
+!!  \param  self  Reference to the dense matrix object
+!!  \param  k     Scalar multiplication factor
+!!  \return new   Reference to the resulting dense matrix object
+!!
+!!  This routine multiplies a dense matrix by a scalar factor.
+! ------------------------------------------------------------------------
+FUNCTION dm_scal_mul_mat_scal_ (self, k) RESULT(new)
+  TYPE(matrixT), INTENT(IN) :: self
+  DOUBLE PRECISION, INTENT(IN) :: k
+  TYPE(matrixT) :: new
+
+  !> initialize solution matrix
+  CALL dm_init(new, dm_numRows(self),dm_numCols(self))
+
+  !> perform multiplication
+  new%dat = k * self%dat
+
+END FUNCTION dm_scal_mul_mat_scal_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Scalar multiplication (k*self version)
+!!
+!!  \param  self  Reference to the dense matrix object
+!!  \param  k     Scalar multiplication factor
+!!  \return new   Reference to the resulting dense matrix object
+!!
+!!  This routine multiplies a dense matrix by a scalar factor. The
+!!  overload is necessary for the OPERATOR (*) overload. This version
+!!  simply calls the other version for better maintainability.
+! ------------------------------------------------------------------------
+FUNCTION dm_scal_mul_scal_mat_ (k, self) RESULT(new)
+  DOUBLE PRECISION, INTENT(IN) :: k
+  TYPE(matrixT), INTENT(IN) :: self
+  TYPE(matrixT) :: new
+
+  new = dm_scalMul(self,k)
+
+END FUNCTION dm_scal_mul_scal_mat_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Post-multiplication of a dense matrix by a vector (non-exception checking)
+!!
+!!  \param  self  Reference to the first (left-hand) dense matrix object
+!!  \param  other Reference to the vector (right-hand) object
+!!  \return new   Reference to the resulting vector object
+!!
+!!  This routine post-multiplies a dense matrix by a vector. That is, it
+!!  performs self*other, not other*self.
+! ------------------------------------------------------------------------
+FUNCTION dm_vec_mul_ (self,other) RESULT(new)
+  TYPE(matrixT), INTENT(IN) :: self
+  TYPE(vectorT), INTENT(IN) :: other
+  TYPE(vectorT) :: new
+
+  !> initialize the solution vector
+  CALL vec_init(new, dm_numRows(self))
+
+  !> perform multiplication (using Fortran built-in for efficiency)
+  new%dat = MATMUL(self%dat,other%dat)
+
+END FUNCTION dm_vec_mul_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Post-multiplication of a dense matrix by a vector (exception checking)
+!!
+!!  \param  self  Reference to the first (left-hand) dense matrix object
+!!  \param  other Reference to the vector (right-hand) object
+!!  \param  exc   Error code
+!!  \return new   Reference to the resulting vector object
+!!
+!!  This routine post-multiplies a dense matrix by a vector. That is, it
+!!  performs self*other, not other*self. If the dimensions of the two
+!!  objects are not compatible (number of columns in self does not equal
+!!  number of rows in other), it returns a DIMEN exception.
+! ------------------------------------------------------------------------
+FUNCTION dm_vec_mul_exc_ (self,other, exc) RESULT(new)
+  TYPE(matrixT), INTENT(IN) :: self
+  TYPE(vectorT), INTENT(IN) :: other
+  INTEGER, INTENT(OUT) :: exc
+  TYPE(vectorT) :: new
+
+  !> make sure dimensions are compatible
+  IF (dm_numCols(self).NE.vec_length(other)) THEN
+    exc=DIMEN
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  new = dm_vecMul(self,other)
+
+END FUNCTION dm_vec_mul_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Multiplication of two dense matrices (non-exception checking)
+!!
+!!  \param  self  Reference to the first (left-hand) dense matrix object
+!!  \param  other Reference to the second (right-hand) dense matrix object
+!!  \return new   Reference to the resulting matrix object
+!!
+!!  This routine multiplies two dense matrices.
+! ------------------------------------------------------------------------
+FUNCTION dm_mat_mul_ (self,other) RESULT(new)
+  TYPE(matrixT), INTENT(IN) :: self, other
+  TYPE(matrixT) :: new
+  INTEGER :: m,n
+
+  !> initialize solution matrix
+  CALL dm_init(new, dm_numRows(self),dm_numCols(other))
+
+  !> perform multiplication (using Fortran built-in for efficiency)
+  new%dat = MATMUL(self%dat,other%dat)
+
+END FUNCTION dm_mat_mul_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Multiplication of two dense matrices (non-exception checking)
+!!
+!!  \param  self  Reference to the first (left-hand) dense matrix object
+!!  \param  other Reference to the second (right-hand) dense matrix object
+!!  \param  exc   Error code
+!!  \return new   Reference to the resulting matrix object
+!!
+!!  This routine multiplies two dense matrices. If the dimensions of the
+!!  two objects are not compatible (number of columns in self does not
+!!  equal number of rows in other), it returns a DIMEN exception.
+! ------------------------------------------------------------------------
+FUNCTION dm_mat_mul_exc_ (self,other, exc) RESULT(new)
+  TYPE(matrixT), INTENT(IN) :: self, other
+  INTEGER, INTENT(OUT) :: exc
+  TYPE(matrixT) :: new
+
+  !> make sure dimensions are compatible
+  IF (dm_numCols(self).NE.dm_numRows(other)) THEN
+    exc=DIMEN
+    CALL log_printLogMsg(exc,sdr)
+    RETURN
+  ELSE
+    exc=OK
+  END IF
+
+  !> call non-exception version
+  new = dm_matMul(self,other)
+
+END FUNCTION dm_mat_mul_exc_
+
+
+! ------------------------------------------------------------------------
+!>  \brief Transpose of a matrix
+!!
+!!  \param  self  Reference to the input dense matrix object
+!!  \return new   Reference to the transposed matrix object
+! ------------------------------------------------------------------------
+FUNCTION dm_transpose_ (self) RESULT (new)
+  TYPE(matrixT), INTENT(IN) :: self
+  TYPE(matrixT) :: new
+
+  !> initialize result matrix
+  CALL dm_init(new, dm_numCols(self),dm_numRows(self))
+
+  !> take the transpose (using Fortran built-in for efficiency)
+  new%dat = TRANSPOSE(self%dat)
+
+END FUNCTION dm_transpose_
+
+END MODULE dense_matrix_def
+
+
+
+! ************************************************************************
+! LOCAL FUNCTIONS
+! ************************************************************************
+
+! none

Examples/SoilWaterStructureInteraction/dynsws/trunk/.svn/text-base/dense_matrix_test.f90.svn-base

+! ------------------------------------------------------------------------
+!>  \brief Module for testing Dense Matrix data type
+! ------------------------------------------------------------------------
+MODULE dense_matrix_test
+USE fruit                   !> Unit testing framework
+USE log_message_control     !> Printing log/error messages
+USE log_messages            !> Log/error codes
+USE dense_matrix_def        !> Dense Matrix data type
+USE vector_def              !> Vector data type
+
+IMPLICIT NONE
+
+
+CONTAINS
+
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for OK exception message on allocation
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test dense matrix
+!!  \param testName Filename for log file (required for exceptions)
+!!  \param m        Number of rows in the matrix
+!!  \param n        Number of columns in the matrix
+!!  \param expMsg   Expected log message
+!!  \param actMsg   Actual log message
+!!
+!!  This test initializes a matrixT and makes sure that the exception
+!!  message is OK (i.e. allocation did not fail).
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_allocation_MSG
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_allocation_MSG'
+  TYPE(matrixT) :: test
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: m=3,n=2
+  INTEGER, PARAMETER :: expMsg=OK
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log message file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize matrix and check the exception
+  CALL dm_init(test, m,n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrix
+  CALL log_closeLogFile()
+  CALL dm_clean(test)
+
+END SUBROUTINE test_dm_allocation_MSG
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for SZE exception message on allocation
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test dense matrix
+!!  \param testName Filename for log file (required for exceptions)
+!!  \param m        Number of rows in the matrix
+!!  \param n        Number of columns in the matrix
+!!  \param expMsg   Expected log message
+!!  \param actMsg   Actual log message
+!!
+!!  This test attempts to initialize a matrixT with invalid size
+!!  parameters and verifies that the correct exception is returned.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_allocation_SZE
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_allocation_SZE'
+  TYPE(matrixT) :: test
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: m=3,n=2
+  INTEGER, PARAMETER :: expMsg=SZE
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> try to initialize with m=-1
+  CALL dm_init(test, -1,n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to initialize with n=-1
+  CALL dm_init(test, m,-1, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrix
+  CALL log_closeLogFile()
+  CALL dm_clean(test)
+
+END SUBROUTINE test_dm_allocation_SZE
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for initialization of matrix data to zeros.
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test dense matrix
+!!  \param m        Number of rows in the matrix
+!!  \param n        Number of columns in the matrix
+!!  \param expDat   Expected initial contents of matrix data
+!!
+!!  This test initializes a matrixT and ensures that the data is
+!!  initialized to zeros.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_allocation_DAT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_allocation_DAT'
+  TYPE(matrixT) :: test
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m,n) :: expDat
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up the expected data
+  expDat = RESHAPE(  (/ 0.d0,0.d0,0.d0,0.d0,0.d0,0.d0 /), SHAPE(expDat) )
+
+  !> initialize matrix and test actual data
+  CALL dm_init(test, m,n)
+  CALL assertEquals(expDat,test%dat, m,n)
+
+  !> deallocate matrix
+  CALL dm_clean(test)
+
+END SUBROUTINE test_dm_allocation_DAT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for deallocation
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test       Test dense matrix
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of rows/columns in the matrix
+!!  \param expBefore  Expected allocation status before deallocation
+!!  \param expAfter   Expected allocation status after deallocation
+!!  \param actBefore  Actual allocation status before deallocation
+!!  \param actAfter   Actual allocation status after deallocation
+!!
+!!  This test initializes a matrixT and ensures that it is
+!!  allocated. It then deallocates the object and ensures that it has
+!!  been deallocated.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_deallocation
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_deallocation'
+  TYPE(matrixT) :: test
+  INTEGER, PARAMETER :: m=3,n=2
+  LOGICAL, PARAMETER :: expBefore = .TRUE.
+  LOGICAL, PARAMETER :: expAfter = .FALSE.
+  LOGICAL :: actBefore, actAfter
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the matrix
+  CALL dm_init(test, m,n)
+
+  !> make sure dat is allocated
+  actBefore = ALLOCATED(test%dat)
+  CALL assertEquals(expBefore,actBefore)
+
+  !> deallocate the matrix
+  CALL dm_clean(test)
+
+  !> make sure dat is deallocated
+  actAfter = ALLOCATED(test%dat)
+  CALL assertEquals(expAfter,actAfter)
+
+END SUBROUTINE test_dm_deallocation
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for number of rows when matrix is not allocated
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test       Test dense matrix
+!!  \param expected   Expected number of rows
+!!  \param actual     Actual number of rows
+!!
+!!  This test makes sure that the number of rows is returned as 0 when
+!!  the matrix is not allocated.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_num_rows_not_allocated
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_num_rows_not_allocated'
+  TYPE(matrixT) :: test
+  INTEGER, PARAMETER :: expected = 0
+  INTEGER :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> check number of rows
+  actual = dm_numRows(test)
+  CALL assertEquals(expected, actual)
+
+END SUBROUTINE test_dm_num_rows_not_allocated
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for number of rows when matrix is allocated
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test       Test dense matrix
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of columns in the matrix
+!!  \param expected   Expected number of rows
+!!  \param actual     Actual number of rows
+!!
+!!  This test checks that the number of rows returned is correct when the
+!!  matrix is allocated.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_num_rows_allocated
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_num_rows_allocated'
+  TYPE(matrixT) :: test
+  INTEGER, PARAMETER :: m=3,n=2
+  INTEGER, PARAMETER :: expected = m
+  INTEGER :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the matrix
+  CALL dm_init(test, m,n)
+
+  !> check number of rows
+  actual = dm_numRows(test)
+  CALL assertEquals(expected, actual)
+
+  !> deallocate the matrix
+  CALL dm_clean(test)
+
+END SUBROUTINE test_dm_num_rows_allocated
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for number of columns when matrix is not allocated
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test       Test dense matrix
+!!  \param expected   Expected number of columns
+!!  \param actual     Actual number of columns
+!!
+!!  This test makes sure that the number of columns is returned as 0 when
+!!  the matrix is not allocated
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_num_cols_not_allocated
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_num_cols_not_allocated'
+  TYPE(matrixT) :: test
+  INTEGER, PARAMETER :: expected = 0
+  INTEGER :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> check number of rows
+  actual = dm_numCols(test)
+  CALL assertEquals(expected, actual)
+
+END SUBROUTINE test_dm_num_cols_not_allocated
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for number of columns when matrix is allocated
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test       Test dense matrix
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of columns in the matrix
+!!  \param expected   Expected number of columns
+!!  \param actual     Actual number of columns
+!!
+!!  This test checks that the number of columns returned is correct when
+!!  the matrix is allocated.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_num_cols_allocated
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_num_cols_allocated'
+  TYPE(matrixT) :: test
+  INTEGER, PARAMETER :: m=3,n=2
+  INTEGER, PARAMETER :: expected = n
+  INTEGER :: actual
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the matrix
+  CALL dm_init(test, m,n)
+
+  !> check number of columns
+  actual = dm_numCols(test)
+  CALL assertEquals(expected, actual)
+
+  !> deallocate the matrix
+  CALL dm_clean(test)
+
+END SUBROUTINE test_dm_num_cols_allocated
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for POSIT exception from dm_get
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test dense matrix
+!!  \param testName Filename for log file (required for exceptions)
+!!  \param m        Number of rows in the matrix
+!!  \param n        Number of columns in the matrix
+!!  \param expMsg   Expected log message
+!!  \param actMsg   Actual log message
+!!  \param v        Dummy variable for get function return
+!!
+!!  This test checks that a POSIT exception is returned when the requested
+!!  location is not inside the matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_get_POSIT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_get_POSIT'
+  TYPE(matrixT) :: test
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: m=3,n=2
+  INTEGER, PARAMETER :: expMsg = POSIT
+  INTEGER :: actMsg
+  DOUBLE PRECISION :: v
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize matrix
+  CALL dm_init(test, m,n)
+
+  !> try to get beyond last row
+  v = dm_get(test, m+1,n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get before first row
+  v = dm_get(test, 0,n, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get after last column
+  v = dm_get(test, m,n+1, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to get before first column
+  v = dm_get(test, m,0, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize the log file and deallocate the matrix
+  CALL log_closeLogFile()
+  CALL dm_clean(test)
+
+END SUBROUTINE test_dm_get_POSIT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct value return from bsm_get
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test dense matrix
+!!  \param m        Number of rows in the matrix
+!!  \param n        Number of rows/columns in the matrix
+!!  \param i        Row index for test location
+!!  \param j        Column index for test location
+!!  \param expVal   Expected data value
+!!  \param actVal   Actual data value
+!!
+!!  This test checks that the correct value is returned from the get
+!!  access program.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_get_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_get_VAL'
+  TYPE(matrixT) :: test
+  INTEGER, PARAMETER :: m=3,n=2
+  INTEGER, PARAMETER :: i=1,j=2
+  DOUBLE PRECISION, PARAMETER :: expVal = 3.d0
+  DOUBLE PRECISION :: actVal
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the matrix
+  CALL dm_init(test, m,n)
+
+  !> set the value manually in the internal data structure
+  test%dat(i,j) = expVal
+
+  !> get the value using the access program
+  actVal = dm_get(test, i,j)
+  CALL assertEquals(expVal,actVal)
+
+  !> deallocate the matrix
+  CALL dm_clean(test)
+
+END SUBROUTINE test_dm_get_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for POSIT exception from dm_set
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test dense matrix
+!!  \param testName Filename for log file (required for exceptions)
+!!  \param m        Number of rows in the matrix
+!!  \param n        Number of columns in the matrix
+!!  \param expMsg   Expected log message
+!!  \param actMsg   Actual log message
+!!  \param v        Dummy variable for set function input
+!!
+!!  This test checks that a POSIT exception is returned when the requested
+!!  location is not inside the matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_set_POSIT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_set_POSIT'
+  TYPE(matrixT) :: test
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: m=3,n=2
+  INTEGER, PARAMETER :: expMsg = POSIT
+  DOUBLE PRECISION, PARAMETER :: v=0.d0
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize the matrix
+  CALL dm_init(test, m,n)
+
+  !> try to set after last row
+  CALL dm_set(test, m+1,n,v, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set before first row
+  CALL dm_set(test, 0,n,v, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set after last column
+  CALL dm_set(test, m,n+1,v, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> try to set before first column
+  CALL dm_set(test, m,0,v, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrix
+  CALL log_closeLogFile()
+  CALL dm_clean(test)
+
+END SUBROUTINE test_dm_set_POSIT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct value return after dm_set
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test     Test dense matrix
+!!  \param m        Number of rows in the matrix
+!!  \param n        Number of columns in the matrix
+!!  \param i        Row index for test location
+!!  \param j        Column index for test location
+!!  \param expVal   Expected data value
+!!  \param actVal   Actual data value
+!!
+!!  This test checks that the correct value is entered by the set
+!!  access program.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_set_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_set_VAL'
+  TYPE(matrixT) :: test
+  INTEGER, PARAMETER :: m=3,n=2
+  INTEGER, PARAMETER :: i=1,j=2
+  DOUBLE PRECISION, PARAMETER :: expVal = 3.d0
+  DOUBLE PRECISION :: actVal
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize the matrix
+  CALL dm_init(test, m,n)
+
+  !> set the value
+  CALL dm_set(test, i,j,expVal)
+  actVal = dm_get(test, i,j)
+  CALL assertEquals(expVal,actVal)
+
+  !> deallocate matrix
+  CALL dm_clean(test)
+
+END SUBROUTINE test_dm_set_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for DIMEN exception when adding matrices
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      First input dense matrix
+!!  \param test2      Second input dense matrix
+!!  \param test3      Output dense matrix
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of columns in the matrix
+!!  \param expMsg     Expected exception message
+!!  \param actMsg     Actual exception message
+!!
+!!  This test checks that a DIMEN exception is raised when matrices with
+!!  incompatible dimensions are added.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_add_DIMEN
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_add_DIMEN'
+  TYPE(matrixT) :: test1, test2, test3
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: m=3,n=2
+  INTEGER, PARAMETER :: expMsg=DIMEN
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize input matrices (data contents not important)
+  CALL dm_init(test1, m,n)
+  CALL dm_init(test2, m+1,n)  !> note different number of rows
+
+  !> attempt to add the matrices
+  test3 = dm_add(test1,test2, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> re-initialize input matrices (data contents not important)
+  CALL dm_init(test1, m,n)
+  CALL dm_init(test2, m,n+1)  !> note different number of columns
+
+  !> attempt to add the matrices
+  test3 = dm_add(test1,test2, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrices
+  CALL log_closeLogFile()
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+  CALL dm_clean(test3)
+
+END SUBROUTINE test_dm_add_DIMEN
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct addition of matrices
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      First input dense matrix
+!!  \param test2      Second input dense matrix
+!!  \param test3      Output dense matrix
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of columns in the matrix
+!!  \param expDat     Expected result of addition
+!!
+!!  This test checks that the correct result is obtained when adding two
+!!  dense matrices.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_add_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_add_VAL'
+  TYPE(matrixT) :: test1, test2, test3
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize expected result
+  expDat = RESHAPE(  (/ 3.d0,3.d0,3.d0,3.d0,3.d0,3.d0 /), SHAPE(expDat) )
+
+  !> initialize first input matrix
+  !! (loop sets matrix to all ones)
+  CALL dm_init(test1, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> initialize second input matrix
+  !! (loop sets matrix to all twos)
+  CALL dm_init(test2, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test2, i,j, 2.d0)
+    END DO
+  END DO
+
+  !> add the matrices
+  test3 = dm_add(test1,test2)
+  CALL assertEquals(expDat,test3%dat, m,n)
+
+  !> deallocate matrices
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+  CALL dm_clean(test3)
+
+END SUBROUTINE test_dm_add_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct addition of matrices ( with OPERATOR (+) )
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      First input dense matrix
+!!  \param test2      Second input dense matrix
+!!  \param test3      Output dense matrix
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of columns in the matrix
+!!  \param expDat     Expected result of addition
+!!
+!!  This test checks that the correct result is obtained when adding two
+!!  dense matrices using the (+) operator.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_add_OP
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_add_OP'
+  TYPE(matrixT) :: test1, test2, test3
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize expected result
+  expDat = RESHAPE(  (/ 3.d0,3.d0,3.d0,3.d0,3.d0,3.d0 /), SHAPE(expDat) )
+
+  !> initialize first input matrix
+  !! (loop sets matrix to all ones)
+  CALL dm_init(test1, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> initialize second input matrix
+  !! (loop sets matrix to all twos)
+  CALL dm_init(test2, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test2, i,j, 2.d0)
+    END DO
+  END DO
+
+  !> add the matrices
+  test3 = test1 + test2
+  CALL assertEquals(expDat,test3%dat, m,n)
+
+  !> deallocate matrices
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+  CALL dm_clean(test3)
+
+END SUBROUTINE test_dm_add_OP
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for scalar multiplication
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1    Input dense matrix
+!!  \param test2    Output dense matrix
+!!  \param k        Scalar coefficient
+!!  \param m        Number of rows in the matrix
+!!  \param n        Number of columns in the matrix
+!!  \param expDat   Expected data contents
+!!
+!!  This test checks that scalar multiplication is performed correctly.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_scal_mul_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_scal_mul_VAL'
+  TYPE(matrixT) :: test1, test2
+  DOUBLE PRECISION, PARAMETER :: k=3.d0
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize expected result
+  expDat = RESHAPE(  (/ 3.d0,3.d0,3.d0,3.d0,3.d0,3.d0 /), SHAPE(expDat) )
+
+  !> initialize matrix
+  !! (loop sets matrix to all ones)
+  CALL dm_init(test1, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> multiply by scalar
+  test2 = dm_scalMul(test1,k)
+  CALL assertEquals(expDat,test2%dat, m,n)
+
+  !> deallocate matrices
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+
+END SUBROUTINE test_dm_scal_mul_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for scalar multiplication (by zero)
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1    Input dense matrix
+!!  \param test2    Output dense matrix
+!!  \param k        Scalar coefficient (zero)
+!!  \param m        Number of rows in the matrix
+!!  \param n        Number of columns in the matrix
+!!  \param expDat   Expected data contents
+!!
+!!  This test checks that scalar multiplication is performed correctly.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_scal_mul_ZERO
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_scal_mul_ZERO'
+  TYPE(matrixT) :: test1, test2
+  DOUBLE PRECISION, PARAMETER :: k=0.d0
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize expected result
+  expDat = RESHAPE(  (/ 0.d0,0.d0,0.d0,0.d0,0.d0,0.d0 /), SHAPE(expDat) )
+
+  !> initialize matrix
+  !! (loop sets matrix to all ones)
+  CALL dm_init(test1, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> multiply by zero
+  test2 = dm_scalMul(test1,k)
+  CALL assertEquals(expDat,test2%dat, m,n)
+
+  !> deallocate matrices
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+
+END SUBROUTINE test_dm_scal_mul_ZERO
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for scalar multiplication ( using OPERATOR (*) )
+!!
+!!  \param unit_name  Name of unit test (for FRUIT)
+!!  \param test1    Input dense matrix
+!!  \param test2    Output dense matrix
+!!  \param k        Scalar coefficient (zero)
+!!  \param m        Number of rows in the matrix
+!!  \param n        Number of columns in the matrix
+!!  \param expDat   Expected data contents
+!!
+!!  This test checks that scalar multiplication is performed correctly
+!!  when using the overloaded (*) operator.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_scal_mul_OP
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_scal_mul_OP'
+  TYPE(matrixT) :: test1, test2
+  DOUBLE PRECISION, PARAMETER :: k=3.d0
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize expected result
+  expDat = RESHAPE(  (/ 3.d0,3.d0,3.d0,3.d0,3.d0,3.d0 /), SHAPE(expDat) )
+
+  !> initialize matrix
+  !! (loop sets matrix to all ones)
+  CALL dm_init(test1, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> multiply by scalar
+  test2 = k*test1
+  CALL assertEquals(expDat,test2%dat, m,n)
+
+  !> clear output matrix
+  CALL dm_clean(test2)
+
+  !> multiply by scalar (switch order of inputs)
+  test2 = test1*k
+  CALL assertEquals(expDat,test2%dat, m,n)
+
+  !> deallocate matrices
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+
+END SUBROUTINE test_dm_scal_mul_OP
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for DIMEN exception when performing matrix-vector multiplication
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      Input dense matrix
+!!  \param test2      Input vector
+!!  \param test3      Dummy output vector
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of columns in the matrix
+!!  \param expMsg     Expected exception message
+!!  \param actMsg     Actual exception message
+!!
+!!  This test checks that a DIMEN exception is raised when matrix-vector
+!!  multiplication is attempted with a vector whose length does not equal
+!!  the number of rows in the matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_vec_mul_DIMEN
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_vec_mul_DIMEN'
+  TYPE(matrixT) :: test1
+  TYPE(vectorT) :: test2, test3
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: m=3,n=2
+  INTEGER, PARAMETER :: expMsg=DIMEN
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> set up matrix and vector (data contents not important)
+  CALL dm_init(test1, m,n)
+  CALL vec_init(test2, m)
+
+  !> attempt matrix-vector multiplication
+  test3 = dm_vecMul(test1,test2, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate objects
+  CALL log_closeLogFile()
+  CALL dm_clean(test1)
+  CALL vec_clean(test2)
+  CALL vec_clean(test3)
+
+END SUBROUTINE test_dm_vec_mul_DIMEN
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct matrix-vector multiplication
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      Input dense matrix
+!!  \param test2      Input vector
+!!  \param test3      Output vector
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of columns in the matrix
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-vector multiplication is performed
+!!  correctly in the general case (input matrix and vector contain data,
+!!  but it is not of a special form).
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_vec_mul_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_vec_mul_VAL'
+  TYPE(matrixT) :: test1
+  TYPE(vectorT) :: test2, test3
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set expected result
+  expDat = RESHAPE(  (/ 9.d0,12.d0,15.d0 /), SHAPE(expDat) )
+
+  !> initialize input matrix
+  !! (loop sets up data)
+  CALL dm_init(test1, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test1, i,j, ( (DBLE(j)-1.d0)*DBLE(m) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> initialize input vector
+  !! (loop sets up data)
+  CALL vec_init(test2, n)
+  DO i = 1,n
+    CALL vec_set(test2, i, DBLE(i))
+  END DO
+
+  !> perform multiplication
+  test3 = dm_vecMul(test1,test2)
+  CALL assertEquals(expDat,test3%dat, m)
+
+  !> deallocate objects
+  CALL dm_clean(test1)
+  CALL vec_clean(test2)
+  CALL vec_clean(test3)
+
+END SUBROUTINE test_dm_vec_mul_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for matrix-vector multiplication (by the zero vector)
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      Input dense matrix
+!!  \param test2      Input vector (zeros)
+!!  \param test3      Output vector
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of columns in the matrix
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-vector multiplication is performed
+!!  correctly when the input vector is all zeros.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_vec_mul_ZERO
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_vec_mul_ZERO'
+  TYPE(matrixT) :: test1
+  TYPE(vectorT) :: test2, test3
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set expected result
+  expDat = RESHAPE(  (/ 0.d0,0.d0,0.d0 /), SHAPE(expDat) )
+
+  !> initialize input matrix
+  !! (loop sets up data)
+  CALL dm_init(test1, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test1, i,j, ( (DBLE(j)-1.d0)*DBLE(m) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> initialize input vector (all zeros)
+  CALL vec_init(test2, n)
+
+  !> perform multiplication
+  test3 = dm_vecMul(test1,test2)
+  CALL assertEquals(expDat,test3%dat, m)
+
+  !> deallocate objects
+  CALL dm_clean(test1)
+  CALL vec_clean(test2)
+  CALL vec_clean(test3)
+
+END SUBROUTINE test_dm_vec_mul_ZERO
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for matrix-vector multiplication (by the identity matrix)
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      Input dense matrix (identity)
+!!  \param test2      Input vector
+!!  \param test3      Output vector
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of columns in the matrix
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-vector multiplication is performed
+!!  correctly when the input matrix is the identity for multiplication.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_vec_mul_IDENT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_vec_mul_IDENT'
+  TYPE(matrixT) :: test1
+  TYPE(vectorT) :: test2, test3
+  INTEGER, PARAMETER :: m=3
+  DOUBLE PRECISION, DIMENSION(m) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize input matrix
+  !! (loop sets up identity matrix)
+  CALL dm_init(test1, m,m)
+  DO j = 1,m
+    DO i = 1,m
+      IF (i.EQ.j) CALL dm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> initialize input vector
+  !! (loop sets data)
+  CALL vec_init(test2, m)
+  DO i = 1,m
+    CALL vec_set(test2, i, DBLE(i))
+  END DO
+
+  !> expected result is the same as the initial vector
+  expDat = test2%dat
+
+  !> perform multiplication
+  test3 = dm_vecMul(test1,test2)
+  CALL assertEquals(expDat,test3%dat, m)
+
+  !> deallocate objects
+  CALL dm_clean(test1)
+  CALL vec_clean(test2)
+  CALL vec_clean(test3)
+
+END SUBROUTINE test_dm_vec_mul_IDENT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct matrix-vector multiplication ( using OPERATOR (*) )
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      Input dense matrix
+!!  \param test2      Input vector
+!!  \param test3      Output vector
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of columns in the matrix
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-vector multiplication is performed
+!!  correctly in the general case (input matrix and vector contain data,
+!!  but it is not of a special form). The overloaded operator (*) is used.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_vec_mul_OP
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_vec_mul_OP'
+  TYPE(matrixT) :: test1
+  TYPE(vectorT) :: test2, test3
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result
+  expDat = RESHAPE(  (/ 9.d0,12.d0,15.d0 /), SHAPE(expDat) )
+
+  !> initialize input matrix
+  !! (loop sets data)
+  CALL dm_init(test1, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test1, i,j, ( (DBLE(j)-1.d0)*DBLE(m) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> initialize input vector
+  !! (loop sets data)
+  CALL vec_init(test2, n)
+  DO i = 1,n
+    CALL vec_set(test2, i, DBLE(i))
+  END DO
+
+  !> perform multiplication
+  test3 = test1*test2
+  CALL assertEquals(expDat,test3%dat, m)
+
+  !> deallocate objects
+  CALL dm_clean(test1)
+  CALL vec_clean(test2)
+  CALL vec_clean(test3)
+
+END SUBROUTINE test_dm_vec_mul_OP
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for DIMEN exception when performing matrix-matrix multiplication
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      First (left-hand) input dense matrix
+!!  \param test2      Second (right-hand) input dense matrix
+!!  \param test3      Dummy output matric
+!!  \param testName   Filename for log file (required for exceptions)
+!!  \param m          Number of rows in the matrix
+!!  \param n          Number of columns in the matrix
+!!  \param expMsg     Expected exception message
+!!  \param actMsg     Actual exception message
+!!
+!!  This test checks that a DIMEN exception is raised when matrix-matrix
+!!  multiplication is attempted between matrices whose dimensions are not
+!!  compatible (i.e. number of columns in test1 does not equal number of
+!!  rows in test2).
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_mat_mul_DIMEN
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_mat_mul_DIMEN'
+  TYPE(matrixT) :: test1, test2, test3
+  CHARACTER (LEN=*), PARAMETER :: testName = 'testName'
+  INTEGER, PARAMETER :: m=3,n=2
+  INTEGER, PARAMETER :: expMsg=DIMEN
+  INTEGER :: actMsg
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize log file
+  CALL log_setFileName(testName)
+  CALL log_initLogFile()
+
+  !> initialize matrices (data content not important)
+  CALL dm_init(test1, m,n)
+  CALL dm_init(test2, m,n)
+
+  !> attempt matrix-matrix multiplication
+  test3 = dm_matMul(test1,test2, actMsg)
+  CALL assertEquals(expMsg,actMsg)
+
+  !> finalize log file and deallocate matrices
+  CALL log_closeLogFile()
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+  CALL dm_clean(test3)
+
+END SUBROUTINE test_dm_mat_mul_DIMEN
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct matrix-matrix multiplication
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      First (left-hand) input dense matrix
+!!  \param test2      Second (right-hand) input dense matrix
+!!  \param test3      Output matrix
+!!  \param m          Number of rows in matrices
+!!  \param n          Number of columns in matrices
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-matrix multiplication is performed
+!!  correctly in the general case (input matrices contain data,
+!!  but it is not of a special form).
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_mat_mul_VAL
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_mat_mul_VAL'
+  TYPE(matrixT) :: test1, test2, test3
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result
+  expDat = RESHAPE( (/ &
+                    30.d0,36.d0,42.d0,66.d0,81.d0,96.d0 &
+                    /), SHAPE(expDat) )
+
+  !> initialize first matrix (dimensions of m x m)
+  !! (loops sets up data)
+  CALL dm_init(test1, m,m)
+  DO j = 1,m
+    DO i = 1,m
+      CALL dm_set(test1, i,j, ( (DBLE(j)-1.d0)*DBLE(m) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> initialize second matrix (dimensions of m x n)
+  !! (loops sets up data)
+  CALL dm_init(test2, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test2, i,j, ( (DBLE(j)-1.d0)*DBLE(m) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> perform matrix-matrix multiplication
+  test3 = dm_matMul(test1,test2)
+  CALL assertEquals(expDat,test3%dat, m,n)
+
+  !> deallocate matrices
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+  CALL dm_clean(test3)
+
+END SUBROUTINE test_dm_mat_mul_VAL
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct matrix-matrix multiplication (by zero matrix)
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      First (left-hand) input dense matrix (zeros)
+!!  \param test2      Second (right-hand) input dense matrix
+!!  \param test3      Output matrix
+!!  \param m          Number of rows in matrices
+!!  \param n          Number of columns in matrices
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-matrix multiplication is performed
+!!  correctly for multiplication by a zero matrix.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_mat_mul_ZERO
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_mat_mul_ZERO'
+  TYPE(matrixT) :: test1, test2, test3
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result
+  expDat = RESHAPE(  (/ 0.d0,0.d0,0.d0,0.d0,0.d0,0.d0 /), SHAPE(expDat) )
+
+  !> initialize first matrix (zeros)
+  CALL dm_init(test1, m,m)
+
+  !> initialize second matrix
+  !! (loop sets up data)
+  CALL dm_init(test2, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test2, i,j, ( (DBLE(j)-1.d0)*DBLE(m) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> perform matrix-matrix multiplication
+  test3 = dm_matMul(test1,test2)
+  CALL assertEquals(expDat,test3%dat, m,n)
+
+  !> deallocate matrices
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+  CALL dm_clean(test3)
+
+END SUBROUTINE test_dm_mat_mul_ZERO
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct matrix-matrix multiplication (by identity matrix)
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      First (left-hand) input dense matrix (identity)
+!!  \param test2      Second (right-hand) input dense matrix
+!!  \param test3      Output matrix
+!!  \param m          Number of rows in matrices
+!!  \param n          Number of columns in matrices
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-matrix multiplication is performed
+!!  correctly for multiplication by the identity matrix for multiplication.
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_mat_mul_IDENT
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_mat_mul_IDENT'
+  TYPE(matrixT) :: test1, test2, test3
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize first matrix
+  !! (loop sets up identity matrix)
+  CALL dm_init(test1, m,m)
+  DO j = 1,m
+    DO i = 1,m
+      IF (i.EQ.j) CALL dm_set(test1, i,j, 1.d0)
+    END DO
+  END DO
+
+  !> initialize second matrix
+  !! (loop sets up data)
+  CALL dm_init(test2, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test2, i,j, ( (DBLE(j)-1.d0)*DBLE(m) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> expected result is the same as the second matrix
+  expDat = test2%dat
+
+  !> perform matrix-matrix multiplication
+  test3 = dm_matMul(test1,test2)
+  CALL assertEquals(expDat,test3%dat, m,n)
+
+  !> deallocate matrices
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+  CALL dm_clean(test3)
+
+END SUBROUTINE test_dm_mat_mul_IDENT
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for correct matrix-matrix multiplication ( using OPERATOR (*) )
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      First (left-hand) input dense matrix
+!!  \param test2      Second (right-hand) input dense matrix
+!!  \param test3      Output matrix
+!!  \param m          Number of rows in matrices
+!!  \param n          Number of columns in matrices
+!!  \param expDat     Expected result
+!!
+!!  This test checks that matrix-matrix multiplication is performed
+!!  correctly in the general case (input matrices contain data,
+!!  but it is not of a special form). This test uses the overloaded
+!!  operator (*).
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_mat_mul_OP
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_mat_mul_OP'
+  TYPE(matrixT) :: test1, test2, test3
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(m,n) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> set up expected result
+  expDat = RESHAPE( (/ &
+                    30.d0,36.d0,42.d0,66.d0,81.d0,96.d0 &
+                    /), SHAPE(expDat) )
+
+  !> initialize first matrix (dimensions of m x m)
+  !! (loop sets up data)
+  CALL dm_init(test1, m,m)
+  DO j = 1,m
+    DO i = 1,m
+      CALL dm_set(test1, i,j, ( (DBLE(j)-1.d0)*DBLE(m) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> initialize second matrix (dimensions of m x n)
+  !! (loop sets up data)
+  CALL dm_init(test2, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test2, i,j, ( (DBLE(j)-1.d0)*DBLE(m) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> perform matrix-matrix multiplication
+  test3 = test1*test2
+  CALL assertEquals(expDat,test3%dat, m,n)
+
+  !> deallocate matrices
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+  CALL dm_clean(test3)
+
+END SUBROUTINE test_dm_mat_mul_OP
+
+
+! ------------------------------------------------------------------------
+!>  \test Test for matrix transpose operation
+!!
+!!  \param unit_name    Name of unit test (for FRUIT)
+!!  \param test1      Input dense matrix
+!!  \param test2      Output dense matrix
+!!  \param m          Number of rows in matrices
+!!  \param n          Number of columns in matrices
+!!  \param expDat     Expected result
+! ------------------------------------------------------------------------
+SUBROUTINE test_dm_transpose
+  CHARACTER (LEN=*), PARAMETER :: unit_name = 'test_dm_transpose'
+  TYPE(matrixT) :: test1, test2
+  INTEGER, PARAMETER :: m=3,n=2
+  DOUBLE PRECISION, DIMENSION(n,m) :: expDat
+  INTEGER :: i,j    !> loop variables
+
+  !> initialize unit test
+  CALL set_unit_name(unit_name)
+
+  !> initialize input matrix
+  !! (loop sets up data)
+  CALL dm_init(test1, m,n)
+  DO j = 1,n
+    DO i = 1,m
+      CALL dm_set(test1, i,j, ( (DBLE(j)-1.d0)*DBLE(m) + DBLE(i) ) )
+    END DO
+  END DO
+
+  !> calculate expected result (using Fortran built-in)
+  expDat = TRANSPOSE(test1%dat)
+
+  !> perform transpose operation
+  test2 = dm_transpose(test1)
+  CALL assertEquals(expDat,test2%dat, n,m)
+
+  !> deallocate matrices
+  CALL dm_clean(test1)
+  CALL dm_clean(test2)
+
+END SUBROUTINE test_dm_transpose
+
+END MODULE dense_matrix_test