stree_equation_routines.f90

Go to the documentation of this file.

MODULE stree_equation_routines

  USE base_routines
  USE basis_routines
  USE boundary_conditions_routines
  USE constants
  USE control_loop_routines
  USE cmiss_mpi  
  USE comp_environment
  USE coordinate_routines
  USE distributed_matrix_vector
  USE domain_mappings
  USE equations_routines
  USE equations_mapping_routines
  USE equations_matrices_routines
  USE equations_set_constants
  USE field_routines
  USE field_io_routines
  USE fluid_mechanics_io_routines
  USE input_output
  USE iso_varying_string
  USE kinds
  USE maths 
  USE matrix_vector
  USE mesh_routines
  USE node_routines
  USE problem_constants
  USE strings
  USE solver_routines
  USE timer
  USE types

#include "macros.h"  

  IMPLICIT NONE

  PRIVATE

  PUBLIC stree_equationssetsolutionmethodset
  PUBLIC stree_equationssetspecificationset
  PUBLIC stree_equationssetsetup
  PUBLIC stree_finite_element_calculate
  PUBLIC stree_pre_solve

CONTAINS 

!
!================================================================================================================================
!

  SUBROUTINE stree_equationssetsolutionmethodset(equationsSet,solutionMethod,err,error,*)

    !Argument variables
    TYPE(equations_set_type), POINTER :: equationsSet
    INTEGER(INTG), INTENT(IN) :: solutionMethod
    INTEGER(INTG), INTENT(OUT) :: err
    TYPE(varying_string), INTENT(OUT) :: error
    !Local Variables
    TYPE(varying_string) :: localError
    
    enters("Stree_EquationsSetSolutionMethodSet",err,error,*999)
    
    IF(ASSOCIATED(equationsset)) THEN
      IF(.NOT.ALLOCATED(equationsset%specification)) THEN
        CALL flagerror("Equations set specification is not allocated.",err,error,*999)
      ELSE IF(SIZE(equationsset%specification,1)/=3) THEN
        CALL flagerror("Equations set specification must have three entries for a structured tree type equations set.", &
          & err,error,*999)
      END IF
      SELECT CASE(equationsset%specification(3))
      CASE(equations_set_stree1d0d_subtype)                                
        SELECT CASE(solutionmethod)
        CASE(equations_set_fem_solution_method)
          equationsset%SOLUTION_METHOD=equations_set_fem_solution_method
        CASE(equations_set_nodal_solution_method)
          CALL flag_error("Not implemented.",err,error,*999)
        CASE(equations_set_bem_solution_method)
          CALL flag_error("Not implemented.",err,error,*999)
        CASE(equations_set_fd_solution_method)
          CALL flag_error("Not implemented.",err,error,*999)
        CASE(equations_set_fv_solution_method)
          CALL flag_error("Not implemented.",err,error,*999)
        CASE(equations_set_gfem_solution_method)
          CALL flag_error("Not implemented.",err,error,*999)
        CASE(equations_set_gfv_solution_method)
          CALL flag_error("Not implemented.",err,error,*999)
        CASE DEFAULT
          localerror="The specified solution method of "//trim(number_to_vstring(solutionmethod,"*",err,error))// &
            & " is invalid."
          CALL flag_error(localerror,err,error,*999)
        END SELECT
      CASE DEFAULT
        localerror="The third equations set specification of "// &
          & trim(number_to_vstring(equationsset%specification(3),"*",err,error))// &
          & " is not valid for a Stree type of a fluid mechanics equations set."
        CALL flag_error(localerror,err,error,*999)
      END SELECT
    ELSE
      CALL flag_error("Equations set is not associated.",err,error,*999)
    ENDIF
       
    exits("Stree_EquationsSetSolutionMethodSet")
    RETURN
999 errorsexits("Stree_EquationsSetSolutionMethodSet",err,error)
    RETURN 1
    
  END SUBROUTINE stree_equationssetsolutionmethodset

!
!================================================================================================================================
!

  SUBROUTINE stree_equationssetspecificationset(equationsSet,specification,err,error,*)

    !Argument variables
    TYPE(equations_set_type), POINTER :: equationsSet
    INTEGER(INTG), INTENT(IN) :: specification(:)
    INTEGER(INTG), INTENT(OUT) :: err
    TYPE(varying_string), INTENT(OUT) :: error
    !Local Variables
    INTEGER(INTG) :: subtype
    TYPE(varying_string) :: localError

    enters("Stree_EquationsSetSpecificationSet",err,error,*999)

    IF(ASSOCIATED(equationsset)) THEN
      IF(SIZE(specification,1)/=3) THEN
        CALL flagerror("Equations set specification must have three entries for a Stree type equations set.", &
          & err,error,*999)
      END IF
      subtype=specification(3)
      SELECT CASE(subtype)
      CASE(equations_set_stree1d0d_subtype)
        !ok
      CASE DEFAULT
        localerror="The third equations set specification of "//trim(numbertovstring(subtype,"*",err,error))// &
          & " is not valid for a Stree type of a fluid mechanics equations set."
        CALL flagerror(localerror,err,error,*999)
      END SELECT
      !Set full specification
      IF(ALLOCATED(equationsset%specification)) THEN
        CALL flagerror("Equations set specification is already allocated.",err,error,*999)
      ELSE
        ALLOCATE(equationsset%specification(3),stat=err)
        IF(err/=0) CALL flagerror("Could not allocate equations set specification.",err,error,*999)
      END IF
      equationsset%specification(1:3)=[equations_set_classical_field_class,equations_set_stree_equation_type,subtype]
    ELSE
      CALL flagerror("Equations set is not associated.",err,error,*999)
    END IF

    exits("Stree_EquationsSetSpecificationSet")
    RETURN
999 errorsexits("Stree_EquationsSetSpecificationSet",err,error)
    RETURN 1
    
  END SUBROUTINE stree_equationssetspecificationset

!
!================================================================================================================================
!

  SUBROUTINE stree_equationssetsetup(equationsSet,equationsSetSetup,err,error,*)

    !Argument variables
    TYPE(equations_set_type), POINTER :: equationsSet
    TYPE(equations_set_setup_type), INTENT(INOUT) :: equationsSetSetup
    INTEGER(INTG), INTENT(OUT) :: err
    TYPE(varying_string), INTENT(OUT) :: error
    !Local Variables
    TYPE(decomposition_type), POINTER :: geometricDecomposition
    TYPE(equations_type), POINTER :: equations
    TYPE(equations_mapping_type), POINTER :: equationsMapping
    TYPE(equations_matrices_type), POINTER :: equationsMatrices
    TYPE(equations_set_materials_type), POINTER :: equationsMaterials
    TYPE(equations_set_equations_set_field_type), POINTER :: equationsEquationsSetField
    TYPE(field_type), POINTER :: equationsSetField
    INTEGER(INTG) :: I,geometricScalingType,geometricMeshComponent,geometricComponentNumber
    INTEGER(INTG) :: dependentFieldNumberOfVariables,dependentFieldNumberOfComponents
    INTEGER(INTG) :: materialsFieldNumberOfVariables,materialsFieldNumberOfComponents
    TYPE(varying_string) :: localError

    enters("Stree_EquationsSetSetup",err,error,*999)

    NULLIFY(equations)
    NULLIFY(equationsmapping)
    NULLIFY(equationsmatrices)
    NULLIFY(equationsmaterials)
    NULLIFY(geometricdecomposition)

    IF(ASSOCIATED(equationsset)) THEN
      IF(.NOT.ALLOCATED(equationsset%specification)) THEN
        CALL flagerror("Equations set specification is not allocated.",err,error,*999)
      ELSE IF(SIZE(equationsset%specification,1)/=3) THEN
        CALL flagerror("Equations set specification must have three entries for a Stree type equations set.", &
          & err,error,*999)
      END IF
      SELECT CASE(equationsset%specification(3))
      CASE(equations_set_stree1d0d_subtype)
        SELECT CASE(equationssetsetup%SETUP_TYPE)
        !-----------------------------------------------------------------
        ! I n i t i a l   s e t u p
        !-----------------------------------------------------------------
        CASE(equations_set_setup_initial_type)
          SELECT CASE(equationsset%specification(3))
          CASE(equations_set_stree1d0d_subtype)
            SELECT CASE(equationssetsetup%ACTION_TYPE)
            CASE(equations_set_setup_start_action)
              CALL stree_equationssetsolutionmethodset(equationsset, &
                & equations_set_fem_solution_method,err,error,*999)
              equationsset%SOLUTION_METHOD=equations_set_fem_solution_method
              equationsequationssetfield=>equationsset%EQUATIONS_SET_FIELD
              IF(equationsequationssetfield%EQUATIONS_SET_FIELD_AUTO_CREATED) THEN
                !Create the auto created equations set field field for SUPG element metrics
                CALL field_create_start(equationssetsetup%FIELD_USER_NUMBER,equationsset%REGION, &
                  & equationsequationssetfield%EQUATIONS_SET_FIELD_FIELD,err,error,*999)
                equationssetfield=>equationsequationssetfield%EQUATIONS_SET_FIELD_FIELD
                CALL field_label_set(equationssetfield,"Equations Set Field",err,error,*999)
                CALL field_type_set_and_lock(equationssetfield,field_general_type,&
                  & err,error,*999)
                CALL field_number_of_variables_set(equationssetfield,1,err,error,*999)
                CALL field_variable_types_set_and_lock(equationssetfield,[field_u_variable_type],err,error,*999)
                CALL field_variable_label_set(equationssetfield,field_u_variable_type,"Stree",err,error,*999)
                CALL field_data_type_set_and_lock(equationssetfield,field_u_variable_type,field_dp_type,err,error,*999)
                CALL field_number_of_components_set_and_lock(equationssetfield,field_u_variable_type,1,err,error,*999)
              ENDIF
            CASE(equations_set_setup_finish_action)
              IF(equationsset%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_AUTO_CREATED) THEN
                CALL field_create_finish(equationsset%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD,err,error,*999)
                CALL field_component_values_initialise(equationsset%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD, &
                 & field_u_variable_type,field_values_set_type,1,1.0_dp,err,error,*999)
              ENDIF
            CASE DEFAULT
              localerror="The action type of "//trim(number_to_vstring(equationssetsetup%ACTION_TYPE, &
                & "*",err,error))// " for a setup type of "//trim(number_to_vstring(equationssetsetup% &
                & setup_type,"*",err,error))// " is not implemented for a Stree equations set."
              CALL flag_error(localerror,err,error,*999)
            END SELECT
          CASE DEFAULT
            localerror="The third equations set specification of "// &
              & trim(number_to_vstring(equationsset%specification(3),"*",err,error))// &
              & " is invalid for a Stree equations set."
            CALL flag_error(localerror,err,error,*999)
          END SELECT
        !-----------------------------------------------------------------
        ! G e o m e t r i c   f i e l d
        !-----------------------------------------------------------------
        CASE(equations_set_setup_geometry_type)
          SELECT CASE(equationsset%specification(3))
          CASE(equations_set_stree1d0d_subtype)
            SELECT CASE(equationssetsetup%ACTION_TYPE)
            CASE(equations_set_setup_start_action)
              equationsequationssetfield=>equationsset%EQUATIONS_SET_FIELD
              equationssetfield=>equationsequationssetfield%EQUATIONS_SET_FIELD_FIELD
              IF(equationsequationssetfield%EQUATIONS_SET_FIELD_AUTO_CREATED) THEN
                CALL field_mesh_decomposition_get(equationsset%GEOMETRY%GEOMETRIC_FIELD,geometricdecomposition,err,error,*999)
                CALL field_mesh_decomposition_set_and_lock(equationssetfield,geometricdecomposition,err,error,*999)
                CALL field_geometric_field_set_and_lock(equationssetfield,equationsset%GEOMETRY%GEOMETRIC_FIELD,err,error,*999)
                CALL field_component_mesh_component_get(equationsset%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, &
                  & 1,geometriccomponentnumber,err,error,*999)                
                CALL field_component_mesh_component_set_and_lock(equationssetfield,field_u_variable_type, &
                  & 1,geometriccomponentnumber,err,error,*999)
                CALL field_component_interpolation_set_and_lock(equationssetfield,field_u_variable_type, &
                  & 1,field_constant_interpolation,err,error,*999)
                !Default the field scaling to that of the geometric field
                CALL field_scaling_type_get(equationsset%GEOMETRY%GEOMETRIC_FIELD,geometricscalingtype,err,error,*999)
                CALL field_scaling_type_set(equationsset%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD, &
                  & geometricscalingtype,err,error,*999)
              ELSE
                !Do nothing
              ENDIF
            CASE(equations_set_setup_finish_action)
              ! do nothing
            CASE DEFAULT
              localerror="The action type of "//trim(number_to_vstring(equationssetsetup%ACTION_TYPE,"*",err,error))// &
                & " for a setup type of "//trim(number_to_vstring(equationssetsetup%SETUP_TYPE,"*",err,error))// &
                & " is invalid for a Stree equation."
              CALL flag_error(localerror,err,error,*999)
            END SELECT
          CASE DEFAULT
            localerror="The third equations set specification of "// &
              & trim(number_to_vstring(equationsset%specification(3),"*",err,error))// &
              & " is invalid for a Stree equations set."
            CALL flag_error(localerror,err,error,*999)
          END SELECT
        !-----------------------------------------------------------------
        ! D e p e n d e n t   f i e l d
        !-----------------------------------------------------------------
        CASE(equations_set_setup_dependent_type)
          SELECT CASE(equationsset%specification(3))
          CASE(equations_set_stree1d0d_subtype)
            dependentfieldnumberofvariables=2
            dependentfieldnumberofcomponents=1
            SELECT CASE(equationssetsetup%ACTION_TYPE)
            !Set start action
            CASE(equations_set_setup_start_action)
              IF(equationsset%DEPENDENT%DEPENDENT_FIELD_AUTO_CREATED) THEN
                !Create the auto created dependent field
                !start field creation with name 'DEPENDENT_FIELD'
                CALL field_create_start(equationssetsetup%FIELD_USER_NUMBER,equationsset%REGION, &
                  & equationsset%DEPENDENT%DEPENDENT_FIELD,err,error,*999)
                !start creation of a new field
                CALL field_type_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD,field_general_type,err,error,*999)
                !label the field
                CALL field_label_set(equationsset%DEPENDENT%DEPENDENT_FIELD,"Dependent Field",err,error,*999)
                !define new created field to be dependent
                CALL field_dependent_type_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD, &
                  & field_dependent_type,err,error,*999)
                !look for decomposition rule already defined
                CALL field_mesh_decomposition_get(equationsset%GEOMETRY%GEOMETRIC_FIELD,geometricdecomposition, &
                  & err,error,*999)
                !apply decomposition rule found on new created field
                CALL field_mesh_decomposition_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD, &
                  & geometricdecomposition,err,error,*999)
                !point new field to geometric field
                CALL field_geometric_field_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD,equationsset%GEOMETRY% &
                  & geometric_field,err,error,*999)
                !set number of variables
                CALL field_number_of_variables_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD, &
                  & dependentfieldnumberofvariables,err,error,*999)
                CALL field_variable_types_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD,[field_u_variable_type, &
                  & field_deludeln_variable_type],err,error,*999)
                !set dimension
                CALL field_dimension_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD,field_u_variable_type, &
                  & field_vector_dimension_type,err,error,*999)
                CALL field_dimension_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD,field_deludeln_variable_type, &
                  & field_vector_dimension_type,err,error,*999)
                !set data type
                CALL field_data_type_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD,field_u_variable_type, &
                  & field_dp_type,err,error,*999)
                CALL field_data_type_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD,field_deludeln_variable_type, &
                  & field_dp_type,err,error,*999)
                ! number of components for 1
                CALL field_number_of_components_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD, &
                  & field_u_variable_type,dependentfieldnumberofcomponents,err,error,*999)
                CALL field_number_of_components_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD, &
                  & field_deludeln_variable_type,dependentfieldnumberofcomponents,err,error,*999)
                CALL field_component_mesh_component_get(equationsset%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, & 
                  & 1,geometricmeshcomponent,err,error,*999)
                !Default to the geometric interpolation setup for U,dUdN
                DO i=1,dependentfieldnumberofcomponents
                  CALL field_component_mesh_component_set(equationsset%DEPENDENT%DEPENDENT_FIELD, & 
                    & field_u_variable_type,i,geometricmeshcomponent,err,error,*999)
                  CALL field_component_mesh_component_set(equationsset%DEPENDENT%DEPENDENT_FIELD, & 
                    & field_deludeln_variable_type,i,geometricmeshcomponent,err,error,*999)
                ENDDO
                SELECT CASE(equationsset%SOLUTION_METHOD)
                !Specify fem solution method
                CASE(equations_set_fem_solution_method)
                  DO i=1,dependentfieldnumberofcomponents
                    CALL field_component_interpolation_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD, &
                      & field_u_variable_type,i,field_node_based_interpolation,err,error,*999)
                    CALL field_component_interpolation_set_and_lock(equationsset%DEPENDENT%DEPENDENT_FIELD, &
                      & field_deludeln_variable_type,i,field_node_based_interpolation,err,error,*999)
                  ENDDO
                  CALL field_scaling_type_get(equationsset%GEOMETRY%GEOMETRIC_FIELD,geometricscalingtype,err,error,*999)
                  CALL field_scaling_type_set(equationsset%DEPENDENT%DEPENDENT_FIELD,geometricscalingtype,err,error,*999)
                CASE DEFAULT
                  localerror="The solution method of " &
                    & //trim(number_to_vstring(equationsset%SOLUTION_METHOD,"*",err,error))// " is invalid."
                  CALL flag_error(localerror,err,error,*999)
                END SELECT
              ELSE 
                !Check the user specified field
                CALL field_type_check(equationssetsetup%FIELD,field_general_type,err,error,*999)
                CALL field_dependent_type_check(equationssetsetup%FIELD,field_dependent_type,err,error,*999)
                CALL field_number_of_variables_check(equationssetsetup%FIELD,dependentfieldnumberofvariables,err,error,*999)
                CALL field_variable_types_check(equationssetsetup%FIELD,[field_u_variable_type, &
                  & field_deludeln_variable_type],err,error,*999)
                CALL field_dimension_check(equationssetsetup%FIELD,field_u_variable_type, & 
                  & field_vector_dimension_type,err,error,*999)
                CALL field_dimension_check(equationssetsetup%FIELD,field_deludeln_variable_type, & 
                  & field_vector_dimension_type,err,error,*999)
                CALL field_data_type_check(equationssetsetup%FIELD,field_u_variable_type,field_dp_type,err,error,*999)
                CALL field_data_type_check(equationssetsetup%FIELD,field_deludeln_variable_type,field_dp_type,err,error,*999)
                CALL field_number_of_components_check(equationssetsetup%FIELD,field_u_variable_type, &
                  & dependentfieldnumberofcomponents,err,error,*999)
                CALL field_number_of_components_check(equationssetsetup%FIELD,field_deludeln_variable_type, &
                  & dependentfieldnumberofcomponents,err,error,*999)
                SELECT CASE(equationsset%SOLUTION_METHOD)
                CASE(equations_set_fem_solution_method)
                  CALL field_component_interpolation_check(equationssetsetup%FIELD,field_u_variable_type,1, &
                    & field_node_based_interpolation,err,error,*999)
                  CALL field_component_interpolation_check(equationssetsetup%FIELD,field_deludeln_variable_type,1, &
                    & field_node_based_interpolation,err,error,*999)
                CASE DEFAULT
                  localerror="The solution method of "//trim(number_to_vstring(equationsset%SOLUTION_METHOD, &
                    & "*",err,error))//" is invalid."
                  CALL flag_error(localerror,err,error,*999)
                END SELECT
              ENDIF
            !Specify finish action
            CASE(equations_set_setup_finish_action)
              IF(equationsset%DEPENDENT%DEPENDENT_FIELD_AUTO_CREATED) THEN
                CALL field_create_finish(equationsset%DEPENDENT%DEPENDENT_FIELD,err,error,*999)
              ENDIF
            CASE DEFAULT
              localerror="The third equations set specification of "// &
                & trim(number_to_vstring(equationsset%specification(3),"*",err,error))// &
                & " is invalid for a Stree equations set."
              CALL flag_error(localerror,err,error,*999)
            END SELECT
          CASE DEFAULT
            localerror="The third equations set specification of "// &
              & trim(number_to_vstring(equationsset%specification(3),"*",err,error))// &
              & " is invalid for a Stree equations set."
            CALL flag_error(localerror,err,error,*999)
          END SELECT
        !-----------------------------------------------------------------
        ! M a t e r i a l s   f i e l d 
        !-----------------------------------------------------------------
        CASE(equations_set_setup_materials_type)
          SELECT CASE(equationsset%specification(3))
          CASE(equations_set_stree1d0d_subtype)
            materialsfieldnumberofvariables=2
            materialsfieldnumberofcomponents=27
            SELECT CASE(equationssetsetup%ACTION_TYPE)
            !Specify start action
            CASE(equations_set_setup_start_action)
              equationsmaterials=>equationsset%MATERIALS
              IF(ASSOCIATED(equationsmaterials)) THEN
                IF(equationsmaterials%MATERIALS_FIELD_AUTO_CREATED) THEN
                  !Create the auto created materials field
                  !start field creation with name 'MATERIAL_FIELD'
                  CALL field_create_start(equationssetsetup%FIELD_USER_NUMBER,equationsset%REGION, & 
                    & equationsset%MATERIALS%MATERIALS_FIELD,err,error,*999)
                  CALL field_type_set_and_lock(equationsmaterials%MATERIALS_FIELD,field_material_type,err,error,*999)
                  !label the field
                  CALL field_label_set(equationsmaterials%MATERIALS_FIELD,"Materials Field",err,error,*999)
                  CALL field_dependent_type_set_and_lock(equationsmaterials%MATERIALS_FIELD,field_independent_type, &
                    & err,error,*999)
                  CALL field_mesh_decomposition_get(equationsset%GEOMETRY%GEOMETRIC_FIELD,geometricdecomposition, & 
                    & err,error,*999)
                  !apply decomposition rule found on new created field
                  CALL field_mesh_decomposition_set_and_lock(equationsset%MATERIALS%MATERIALS_FIELD, & 
                    & geometricdecomposition,err,error,*999)
                  !point new field to geometric field
                  CALL field_geometric_field_set_and_lock(equationsmaterials%MATERIALS_FIELD,equationsset%GEOMETRY% &
                    & geometric_field,err,error,*999)
                  CALL field_number_of_variables_set(equationsmaterials%MATERIALS_FIELD, & 
                    & materialsfieldnumberofvariables,err,error,*999)
                  CALL field_variable_types_set_and_lock(equationsmaterials%MATERIALS_FIELD, & 
                    & [field_u_variable_type,field_v_variable_type],err,error,*999)
                  CALL field_dimension_set_and_lock(equationsmaterials%MATERIALS_FIELD,field_u_variable_type, &
                    & field_vector_dimension_type,err,error,*999)
                  CALL field_dimension_set_and_lock(equationsmaterials%MATERIALS_FIELD,field_v_variable_type, &
                    & field_vector_dimension_type,err,error,*999)
                  CALL field_data_type_set_and_lock(equationsmaterials%MATERIALS_FIELD,field_u_variable_type, &
                    & field_dp_type,err,error,*999)
                  CALL field_data_type_set_and_lock(equationsmaterials%MATERIALS_FIELD,field_v_variable_type, &
                    & field_dp_type,err,error,*999)
                  CALL field_number_of_components_set_and_lock(equationsmaterials%MATERIALS_FIELD,field_u_variable_type, &
                    & materialsfieldnumberofcomponents,err,error,*999)
                  CALL field_number_of_components_set_and_lock(equationsmaterials%MATERIALS_FIELD,field_v_variable_type, &
                    & materialsfieldnumberofcomponents,err,error,*999)
                  CALL field_component_mesh_component_get(equationsset%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, &
                    & 1,geometriccomponentnumber,err,error,*999)
                  DO i=1,materialsfieldnumberofcomponents
                    CALL field_component_mesh_component_set(equationsmaterials%MATERIALS_FIELD,field_u_variable_type, &
                      & i,geometriccomponentnumber,err,error,*999)
                    CALL field_component_mesh_component_set(equationsmaterials%MATERIALS_FIELD,field_v_variable_type, &
                      & i,geometriccomponentnumber,err,error,*999)
                  ENDDO
                  DO i=1,materialsfieldnumberofcomponents
                    CALL field_component_interpolation_set(equationsmaterials%MATERIALS_FIELD,field_u_variable_type, &
                      & i,field_node_based_interpolation,err,error,*999)
                    CALL field_component_interpolation_set(equationsmaterials%MATERIALS_FIELD,field_v_variable_type, &
                      & i,field_node_based_interpolation,err,error,*999)
                  ENDDO
                  !Default the field scaling to that of the geometric field
                  CALL field_scaling_type_get(equationsset%GEOMETRY%GEOMETRIC_FIELD,geometricscalingtype,err,error,*999)
                  CALL field_scaling_type_set(equationsmaterials%MATERIALS_FIELD,geometricscalingtype,err,error,*999)
                ELSE
                  !Check the user specified field
                  CALL field_type_check(equationssetsetup%FIELD,field_material_type,err,error,*999)
                  CALL field_dependent_type_check(equationssetsetup%FIELD,field_independent_type,err,error,*999)
                  CALL field_number_of_variables_check(equationssetsetup%FIELD,materialsfieldnumberofvariables,err,error,*999)
                  CALL field_variable_types_check(equationssetsetup%FIELD,[field_u_variable_type,field_u_variable_type], &
                    & err,error,*999)
                  CALL field_dimension_check(equationssetsetup%FIELD,field_u_variable_type,field_vector_dimension_type, &
                    & err,error,*999)
                  CALL field_dimension_check(equationssetsetup%FIELD,field_v_variable_type,field_vector_dimension_type, &
                    & err,error,*999)
                  CALL field_data_type_check(equationssetsetup%FIELD,field_u_variable_type,field_dp_type,err,error,*999)
                  CALL field_data_type_check(equationssetsetup%FIELD,field_v_variable_type,field_dp_type,err,error,*999)
                  CALL field_number_of_components_check(equationssetsetup%FIELD,field_u_variable_type, &
                    & materialsfieldnumberofcomponents,err,error,*999)
                  CALL field_number_of_components_check(equationssetsetup%FIELD,field_v_variable_type, &
                    & materialsfieldnumberofcomponents,err,error,*999)
                ENDIF
              ELSE
                CALL flag_error("Equations set materials is not associated.",err,error,*999)
              END IF
              !Specify start action
            CASE(equations_set_setup_finish_action)
              equationsmaterials=>equationsset%MATERIALS
              IF(ASSOCIATED(equationsmaterials)) THEN
                IF(equationsmaterials%MATERIALS_FIELD_AUTO_CREATED) THEN
                  CALL field_create_finish(equationsmaterials%MATERIALS_FIELD,err,error,*999)
                ENDIF
              ELSE
                CALL flag_error("Equations set materials is not associated.",err,error,*999)
              ENDIF
            CASE DEFAULT
              localerror="The action type of "//trim(number_to_vstring(equationssetsetup%ACTION_TYPE,"*", & 
                & err,error))//" for a setup type of "//trim(number_to_vstring(equationssetsetup%SETUP_TYPE,"*", & 
                & err,error))//" is invalid for Stree equation."
              CALL flag_error(localerror,err,error,*999)
            END SELECT
          CASE DEFAULT
            localerror="The third equations set specification of "// &
              & trim(number_to_vstring(equationsset%specification(3),"*",err,error))// &
              & " is invalid for a Stree equation."
            CALL flag_error(localerror,err,error,*999)
          END SELECT
        !-----------------------------------------------------------------
        ! E q u a t i o n s    t y p e
        !-----------------------------------------------------------------
        CASE(equations_set_setup_equations_type)
          SELECT CASE(equationsset%specification(3))
          CASE(equations_set_stree1d0d_subtype)
            SELECT CASE(equationssetsetup%ACTION_TYPE)
            CASE(equations_set_setup_start_action)
              equationsmaterials=>equationsset%MATERIALS
              IF(ASSOCIATED(equationsmaterials)) THEN              
                IF(equationsmaterials%MATERIALS_FINISHED) THEN
                  CALL equations_create_start(equationsset,equations,err,error,*999)
                  CALL equations_linearity_type_set(equations,equations_linear,err,error,*999)
                  CALL equations_time_dependence_type_set(equations,equations_static,err,error,*999)
                ELSE
                  CALL flag_error("Equations set materials has not been finished.",err,error,*999)
                ENDIF
              ELSE
                CALL flag_error("Equations materials is not associated.",err,error,*999)
              ENDIF
            CASE(equations_set_setup_finish_action)
              SELECT CASE(equationsset%SOLUTION_METHOD)
              CASE(equations_set_fem_solution_method)
                !Finish the creation of the equations
                CALL equations_set_equations_get(equationsset,equations,err,error,*999)
                CALL equations_create_finish(equations,err,error,*999)
                !Create the equations mapping.
                CALL equations_mapping_create_start(equations,equationsmapping,err,error,*999)
                CALL equationsmapping_linearmatricesnumberset(equationsmapping,1,err,error,*999)
                CALL equationsmapping_linearmatricesvariabletypesset(equationsmapping,[field_u_variable_type],err,error,*999)
                CALL equations_mapping_rhs_variable_type_set(equationsmapping,field_deludeln_variable_type,err,error,*999)
                CALL equations_mapping_create_finish(equationsmapping,err,error,*999)
                !Create the equations matrices
                CALL equations_matrices_create_start(equations,equationsmatrices,err,error,*999)
                SELECT CASE(equations%SPARSITY_TYPE)
                CASE(equations_matrices_full_matrices)
                  CALL equations_matrices_linear_storage_type_set(equationsmatrices,[matrix_block_storage_type],err,error,*999)
                CASE(equations_matrices_sparse_matrices)
                  CALL equations_matrices_linear_storage_type_set(equationsmatrices, & 
                    & [matrix_compressed_row_storage_type],err,error,*999)
                  CALL equationsmatrices_linearstructuretypeset(equationsmatrices, & 
                    & [equations_matrix_fem_structure],err,error,*999)
                CASE DEFAULT
                  localerror="The equations matrices sparsity type of "// &
                    & trim(number_to_vstring(equations%SPARSITY_TYPE,"*",err,error))//" is invalid."
                  CALL flag_error(localerror,err,error,*999)
                END SELECT
                CALL equations_matrices_create_finish(equationsmatrices,err,error,*999)
              CASE(equations_set_bem_solution_method)
                CALL flag_error("Not implemented.",err,error,*999)
              CASE(equations_set_fd_solution_method)
                CALL flag_error("Not implemented.",err,error,*999)
              CASE(equations_set_fv_solution_method)
                CALL flag_error("Not implemented.",err,error,*999)
              CASE(equations_set_gfem_solution_method)
                CALL flag_error("Not implemented.",err,error,*999)
              CASE(equations_set_gfv_solution_method)
                CALL flag_error("Not implemented.",err,error,*999)
              CASE DEFAULT
                localerror="The solution method of "//trim(number_to_vstring(equationsset%SOLUTION_METHOD, &
                  & "*",err,error))//" is invalid."
                CALL flag_error(localerror,err,error,*999)
              END SELECT
            CASE DEFAULT
              localerror="The action type of "//trim(number_to_vstring(equationssetsetup%ACTION_TYPE,"*",err,error))// &
                & " for a setup type of "//trim(number_to_vstring(equationssetsetup%SETUP_TYPE,"*",err,error))// &
                & " is invalid for a Strees equation."
              CALL flag_error(localerror,err,error,*999)
            END SELECT
          CASE DEFAULT
            localerror="The third equations set specification of "// &
              & trim(number_to_vstring(equationsset%specification(3),"*",err,error))// &
              & " is invalid for a Stree equation."
            CALL flag_error(localerror,err,error,*999)
          END SELECT
        CASE DEFAULT
          localerror="The setup type of "//trim(number_to_vstring(equationssetsetup%SETUP_TYPE,"*",err,error))// &
            & " is invalid for a Strees equation set."
          CALL flag_error(localerror,err,error,*999)
        END SELECT
      CASE DEFAULT
        localerror="The third equations set specification of "// &
          & trim(number_to_vstring(equationsset%specification(3),"*",err,error))// &
          & " does not equal a Strees equation set."
        CALL flag_error(localerror,err,error,*999)
      END SELECT
    ELSE
      CALL flag_error("Equations set is not associated.",err,error,*999)
    ENDIF

    exits("Stree_EquationsSetSetup")
    RETURN
999 errorsexits("Stree_EquationsSetSetup",err,error)
    RETURN 1
    
  END SUBROUTINE stree_equationssetsetup

  !
  !================================================================================================================================
  !

  SUBROUTINE stree_finite_element_calculate(equationsSet,nodeNumber,err,error,*)

    !Argument variables
    TYPE(equations_set_type), POINTER :: equationsSet
    INTEGER(INTG), INTENT(IN) :: nodeNumber
    INTEGER(INTG), INTENT(OUT) :: err
    TYPE(varying_string), INTENT(OUT) :: error
    !Local Variables
    TYPE(domain_nodes_type), POINTER :: domainNodes
    TYPE(domain_type), POINTER :: domain
    TYPE(equations_type), POINTER :: equations
    TYPE(equations_mapping_type), POINTER :: equationsMapping
    TYPE(equations_matrices_type), POINTER :: equationsMatrices
    TYPE(equations_mapping_linear_type), POINTER :: linearMapping
    TYPE(equations_matrices_linear_type), POINTER :: linearMatrices
    TYPE(equations_matrix_type), POINTER :: stiffnessMatrix
    TYPE(field_type), POINTER :: materialsField,dependentField
    TYPE(field_variable_type), POINTER :: fieldVariable
    TYPE(varying_string) :: localError
    REAL(DP), POINTER :: dependentParameters(:),materialsParameters(:)

    enters("STREE_FINITE_ELEMENT_CALCULATE",err,error,*999)

    NULLIFY(equations)
    NULLIFY(equationsmapping)
    NULLIFY(equationsmatrices)
    NULLIFY(linearmapping)
    NULLIFY(linearmatrices)
    NULLIFY(stiffnessmatrix)
    NULLIFY(dependentfield)
    NULLIFY(materialsfield)
    NULLIFY(domain)
    NULLIFY(domainnodes)
    NULLIFY(dependentparameters)
    NULLIFY(materialsparameters)
    NULLIFY(fieldvariable)

    IF(ASSOCIATED(equationsset)) THEN
      equations=>equationsset%EQUATIONS
      IF(ASSOCIATED(equations)) THEN
        dependentfield=>equations%EQUATIONS_SET%DEPENDENT%DEPENDENT_FIELD
        IF(ASSOCIATED(dependentfield)) THEN
          domain=>dependentfield%DECOMPOSITION%DOMAIN(dependentfield%DECOMPOSITION%MESH_COMPONENT_NUMBER)%PTR
          IF(ASSOCIATED(domain)) THEN
            domainnodes=>domain%TOPOLOGY%NODES
          ELSE
            CALL flag_error("Domain is not associated.",err,error,*999)
          ENDIF
        ELSE
          CALL flag_error("Dependent Field is not associated.",err,error,*999)
        ENDIF
      ELSE
        CALL flag_error("Equations set equations is not associated.",err,error,*999)
      ENDIF
    ELSE
      CALL flag_error("Equations set is not associated.",err,error,*999)
    ENDIF

    SELECT CASE(equationsset%specification(3))
    CASE(equations_set_stree1d0d_subtype)
      !Set General and Specific Pointers
      equationsmatrices=>equations%EQUATIONS_MATRICES
      equationsmapping=>equations%EQUATIONS_MAPPING
      linearmatrices=>equationsmatrices%LINEAR_MATRICES
      stiffnessmatrix=>linearmatrices%MATRICES(1)%PTR
      linearmapping=>equationsmapping%LINEAR_MAPPING
      stiffnessmatrix%ELEMENT_MATRIX%matrix=0.0_dp

    CASE DEFAULT
      localerror="The third equations set specification of "// &
        & trim(number_to_vstring(equationsset%specification(3),"*",err,error))// &
        & " is not valid for a Stree type of a fluid mechanics equations set."
      CALL flag_error(localerror,err,error,*999)
    END SELECT

    exits("STREE_FINITE_ELEMENT_CALCULATE")
    RETURN
999 errorsexits("STREE_FINITE_ELEMENT_CALCULATE",err,error)
    RETURN 1
    
  END SUBROUTINE stree_finite_element_calculate

  !
  !================================================================================================================================
  !      

  SUBROUTINE stree_pre_solve(solver,err,error,*)

    !Argument variables
    TYPE(solver_type), POINTER :: solver
    INTEGER(INTG), INTENT(OUT) :: err
    TYPE(varying_string), INTENT(OUT) :: error
    !Local Variables
    TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
    TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
    TYPE(control_loop_type), POINTER :: controlLoop,parentLoop,navierstokesLoop
    TYPE(domain_nodes_type), POINTER :: domainNodes
    TYPE(domain_type), POINTER :: domain
    TYPE(equations_set_type), POINTER :: equationsSet,navierstokesEquationsSet
    TYPE(equations_type), POINTER :: equations,navierstokesEquations
    TYPE(field_type), POINTER :: materialsField,navierstokesDependentField
    TYPE(field_variable_type), POINTER :: dependentFieldVariable
    TYPE(solver_equations_type), POINTER :: solverEquations,navierstokesSolverEquations
    TYPE(solver_mapping_type), POINTER :: solverMapping,navierstokesSolverMapping
    TYPE(solvers_type), POINTER :: solvers
    TYPE(solver_type), POINTER :: navierstokesSolver
    TYPE(varying_string) :: localError
    INTEGER(INTG) :: variableIdx,BOUNDARY_CONDITION_CHECK_VARIABLE,nodeIdx,componentIdx,derivativeIdx,versionIdx
    INTEGER(INTG) :: dependentDof,dependentVariableType,userNodeNumber,m
    REAL(DP) :: currentTime,timeIncrement,flow

    enters("Stree_PRE_SOLVE",err,error,*999)

    ! Some preliminary sanity checks
    IF(ASSOCIATED(solver)) THEN
      solvers=>solver%SOLVERS
      IF(ASSOCIATED(solvers)) THEN
        controlloop=>solvers%CONTROL_LOOP
        CALL control_loop_current_times_get(controlloop,currenttime,timeincrement,err,error,*999)
        parentloop=>controlloop%PARENT_LOOP
        navierstokesloop=>parentloop%SUB_LOOPS(2)%PTR
        navierstokessolver=>navierstokesloop%SOLVERS%SOLVERS(2)%PTR
        IF(ASSOCIATED(controlloop%PROBLEM)) THEN
          SELECT CASE(controlloop%PROBLEM%specification(3))
          CASE(problem_stree1d0d_navier_stokes_subtype, &
             & problem_stree1d0d_adv_navier_stokes_subtype)
            solverequations=>solver%SOLVER_EQUATIONS
            navierstokessolverequations=>navierstokessolver%SOLVER_EQUATIONS
            IF(ASSOCIATED(solverequations)) THEN
              solvermapping=>solverequations%SOLVER_MAPPING
              navierstokessolvermapping=>navierstokessolverequations%SOLVER_MAPPING
              IF(ASSOCIATED(solvermapping)) THEN
                equationsset=>solvermapping%EQUATIONS_SETS(1)%PTR
                navierstokesequationsset=>navierstokessolvermapping%EQUATIONS_SETS(1)%PTR
                IF(ASSOCIATED(equationsset)) THEN
                  equations=>equationsset%EQUATIONS
                  navierstokesequations=>navierstokesequationsset%EQUATIONS
                  IF(ASSOCIATED(equations)) THEN
                    materialsfield=>equationsset%MATERIALS%MATERIALS_FIELD
                    navierstokesdependentfield=>navierstokesequationsset%DEPENDENT%DEPENDENT_FIELD
                    IF(.NOT.ASSOCIATED(materialsfield)) THEN
                      CALL flag_error("Dependent field is not associated.",err,error,*999)
                    END IF
                  ELSE
                    CALL flag_error("Equations set equations is not associated.",err,error,*999)
                  END IF
                ELSE
                  CALL flag_error("Equations set is not associated.",err,error,*999)
                END IF
              ELSE
                CALL flag_error("Solver mapping is not associated.",err,error,*999)
              END IF
            ELSE
              CALL flag_error("Solver equations is not associated.",err,error,*999)
            END IF
          CASE DEFAULT
            localerror="The third problem specification of "// &
              & trim(number_to_vstring(controlloop%PROBLEM%specification(3),"*",err,error))// &
              & " is not valid for boundary flux calculation."
            CALL flag_error(localerror,err,error,*999)
          END SELECT
        ELSE
          CALL flag_error("Problem is not associated.",err,error,*999)
        END IF
      ELSE
        CALL flag_error("Solvers is not associated.",err,error,*999)
      END IF
    ELSE
      CALL flag_error("Solver is not associated.",err,error,*999)
    END IF

    boundary_conditions=>navierstokessolverequations%BOUNDARY_CONDITIONS
    DO variableidx=1,navierstokesdependentfield%NUMBER_OF_VARIABLES
      dependentvariabletype=navierstokesdependentfield%VARIABLES(variableidx)%VARIABLE_TYPE
      NULLIFY(dependentfieldvariable)
      CALL field_variable_get(navierstokesdependentfield,dependentvariabletype,dependentfieldvariable,err,error,*999)
      CALL boundary_conditions_variable_get(boundary_conditions, &
        & dependentfieldvariable,boundary_conditions_variable,err,error,*999)
      IF(ASSOCIATED(boundary_conditions_variable)) THEN
        IF(ASSOCIATED(dependentfieldvariable)) THEN
          DO componentidx=1,dependentfieldvariable%NUMBER_OF_COMPONENTS
            IF(dependentfieldvariable%COMPONENTS(componentidx)%INTERPOLATION_TYPE==field_node_based_interpolation) THEN
              domain=>dependentfieldvariable%COMPONENTS(componentidx)%DOMAIN
              IF(ASSOCIATED(domain)) THEN
                IF(ASSOCIATED(domain%TOPOLOGY)) THEN
                  domainnodes=>domain%TOPOLOGY%NODES
                  IF(ASSOCIATED(domainnodes)) THEN
                    !Loop over the local nodes excluding the ghosts.
                    DO nodeidx=1,domainnodes%NUMBER_OF_NODES
                      usernodenumber=domainnodes%NODES(nodeidx)%USER_NUMBER
                      DO derivativeidx=1,domainnodes%NODES(nodeidx)%NUMBER_OF_DERIVATIVES
                        DO versionidx=1,domainnodes%NODES(nodeidx)%DERIVATIVES(derivativeidx)%numberOfVersions
                          dependentdof = dependentfieldvariable%COMPONENTS(componentidx)%PARAM_TO_DOF_MAP% &
                            & node_param2dof_map%NODES(nodeidx)%DERIVATIVES(derivativeidx)%VERSIONS(versionidx)
                          boundary_condition_check_variable=boundary_conditions_variable%CONDITION_TYPES(dependentdof)
                          IF(boundary_condition_check_variable==boundary_condition_fixed_stree) THEN  
                            ! Update dependent field value
                            IF(ASSOCIATED(materialsfield)) THEN
                              CALL field_parameter_set_get_node(navierstokesdependentfield,field_u_variable_type, &
                                & field_values_set_type,versionidx,derivativeidx,nodeidx,1,flow,err,error,*999)
                              m=int(currenttime)-800*(int(currenttime)/800)
                              CALL field_parameter_set_update_local_node(materialsfield,field_v_variable_type, &
                                & field_values_set_type,versionidx,derivativeidx,m+1,1,flow,err,error,*999)
                            ENDIF
                          ENDIF
                        ENDDO !versionIdx
                      ENDDO !derivativeIdx
                    ENDDO !nodeIdx
                  ELSE
                    CALL flag_error("Domain topology nodes is not associated.",err,error,*999)
                  ENDIF
                ELSE
                  CALL flag_error("Domain topology is not associated.",err,error,*999)
                ENDIF
              ELSE
                CALL flag_error("Domain is not associated.",err,error,*999)
              ENDIF
            ELSE
              CALL flag_error("Only node based interpolation is implemented.",err,error,*999)
            ENDIF
          ENDDO !componentIdx
        ELSE
          CALL flag_error("Dependent field variable is not associated.",err,error,*999)
        ENDIF
      ENDIF
    ENDDO !variableIdx

    exits("Stree_PRE_SOLVE")
    RETURN
999 errorsexits("Stree_PRE_SOLVE",err,error)
    RETURN 1
    
  END SUBROUTINE stree_pre_solve

  !
  !================================================================================================================================
  !      

END MODULE stree_equation_routines