diff --git a/src/CPFEM.f90 b/src/CPFEM.f90
index 713aab5d7..8a12bef76 100644
--- a/src/CPFEM.f90
+++ b/src/CPFEM.f90
@@ -19,7 +19,7 @@ module CPFEM
   use HDF5_utilities
   use results
   use lattice
-  use constitutive
+  use phase
 
   implicit none
   private
diff --git a/src/CPFEM2.f90 b/src/CPFEM2.f90
index e696858cf..2b32a0cbb 100644
--- a/src/CPFEM2.f90
+++ b/src/CPFEM2.f90
@@ -19,7 +19,7 @@ module CPFEM2
   use discretization
   use HDF5_utilities
   use homogenization
-  use constitutive
+  use phase
 #if    defined(Mesh)
   use FEM_quadrature
   use discretization_mesh
diff --git a/src/commercialFEM_fileList.f90 b/src/commercialFEM_fileList.f90
index 311129d61..d145d3965 100644
--- a/src/commercialFEM_fileList.f90
+++ b/src/commercialFEM_fileList.f90
@@ -32,15 +32,15 @@
 #include "phase_mechanics_eigendeformation.f90"
 #include "phase_mechanics_eigendeformation_cleavageopening.f90"
 #include "phase_mechanics_eigendeformation_slipplaneopening.f90"
-#include "phase_thermal.f90"
-#include "phase_mechanics_eigendeformation_thermalexpansion.f90"
-#include "phase_thermal_dissipation.f90"
-#include "phase_thermal_externalheat.f90"
 #include "phase_damage.f90"
 #include "phase_damage_isobrittle.f90"
 #include "phase_damage_isoductile.f90"
 #include "phase_damage_anisobrittle.f90"
 #include "phase_damage_anisoductile.f90"
+#include "phase_thermal.f90"
+#include "phase_mechanics_eigendeformation_thermalexpansion.f90"
+#include "phase_thermal_dissipation.f90"
+#include "phase_thermal_externalheat.f90"
 #include "damage_none.f90"
 #include "damage_nonlocal.f90"
 #include "homogenization.f90"
diff --git a/src/damage_nonlocal.f90 b/src/damage_nonlocal.f90
index f566ebbeb..807231889 100644
--- a/src/damage_nonlocal.f90
+++ b/src/damage_nonlocal.f90
@@ -8,7 +8,7 @@ module damage_nonlocal
   use config
   use YAML_types
   use lattice
-  use constitutive
+  use phase
   use results
 
   implicit none
diff --git a/src/homogenization.f90 b/src/homogenization.f90
index 9b8c33c2e..4738b3ad8 100644
--- a/src/homogenization.f90
+++ b/src/homogenization.f90
@@ -10,7 +10,7 @@ module homogenization
   use config
   use math
   use material
-  use constitutive
+  use phase
   use discretization
   use damage_none
   use damage_nonlocal
diff --git a/src/phase.f90 b/src/phase.f90
index 65fb7a859..849c5d8c6 100644
--- a/src/phase.f90
+++ b/src/phase.f90
@@ -3,7 +3,7 @@
 !> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH
 !> @brief elasticity, plasticity, damage & thermal internal microstructure state
 !--------------------------------------------------------------------------------------------------
-module constitutive
+module phase
   use prec
   use math
   use rotations
@@ -16,7 +16,6 @@ module constitutive
   use parallelization
   use HDF5_utilities
   use results
-
   implicit none
   private
 
@@ -304,14 +303,14 @@ module constitutive
         dLd_dTstar                                                                                  !< derivative of Ld with respect to Tstar (4th-order tensor)
     end subroutine kinematics_slipplane_opening_LiAndItsTangent
 
-    module subroutine kinematics_thermal_expansion_LiAndItsTangent(Li, dLi_dTstar, ph,me)
+    module subroutine thermalexpansion_LiAndItsTangent(Li, dLi_dTstar, ph,me)
       integer, intent(in) :: ph, me
                                                                                        !< element number
       real(pReal),   intent(out), dimension(3,3) :: &
         Li                                                                                          !< thermal velocity gradient
       real(pReal),   intent(out), dimension(3,3,3,3) :: &
         dLi_dTstar                                                                                  !< derivative of Li with respect to Tstar (4th-order tensor defined to be zero)
-    end subroutine kinematics_thermal_expansion_LiAndItsTangent
+    end subroutine thermalexpansion_LiAndItsTangent
 
     module subroutine plastic_dependentState(co,ip,el)
       integer, intent(in) :: &
@@ -325,6 +324,21 @@ module constitutive
 
 
   type(tDebugOptions) :: debugConstitutive
+#if __INTEL_COMPILER >= 1900
+  public :: &
+    prec, &
+    math, &
+    rotations, &
+    IO, &
+    config, &
+    material, &
+    results, &
+    lattice, &
+    discretization, &
+    parallelization, &
+    HDF5_utilities, &
+    results
+#endif
 
   public :: &
     constitutive_init, &
@@ -788,4 +802,4 @@ subroutine constitutive_restartRead(fileHandle)
 end subroutine constitutive_restartRead
 
 
-end module constitutive
+end module phase
diff --git a/src/phase_damage.f90 b/src/phase_damage.f90
index 72be51be1..6c1a0b8d0 100644
--- a/src/phase_damage.f90
+++ b/src/phase_damage.f90
@@ -1,7 +1,7 @@
 !----------------------------------------------------------------------------------------------------
 !> @brief internal microstructure state for all damage sources and kinematics constitutive models
 !----------------------------------------------------------------------------------------------------
-submodule(constitutive) constitutive_damage
+submodule(phase) damagee
   enum, bind(c); enumerator :: &
     DAMAGE_UNDEFINED_ID, &
     DAMAGE_ISOBRITTLE_ID, &
@@ -525,4 +525,4 @@ module function constitutive_damage_get_phi(co,ip,el) result(phi)
 end function constitutive_damage_get_phi
 
 
-end submodule constitutive_damage
+end submodule damagee
diff --git a/src/phase_damage_anisobrittle.f90 b/src/phase_damage_anisobrittle.f90
index 095046c1c..f78a0ec99 100644
--- a/src/phase_damage_anisobrittle.f90
+++ b/src/phase_damage_anisobrittle.f90
@@ -4,7 +4,7 @@
 !> @brief material subroutine incorporating anisotropic brittle damage source mechanism
 !> @details to be done
 !--------------------------------------------------------------------------------------------------
-submodule (constitutive:constitutive_damage) source_damage_anisoBrittle
+submodule (phase:damagee) anisobrittle
 
   integer,                       dimension(:),           allocatable :: &
     source_damage_anisoBrittle_offset, &                                                            !< which source is my current source mechanism?
@@ -215,4 +215,4 @@ module subroutine anisobrittle_results(phase,group)
 
 end subroutine anisobrittle_results
 
-end submodule source_damage_anisoBrittle
+end submodule anisobrittle
diff --git a/src/phase_damage_anisoductile.f90 b/src/phase_damage_anisoductile.f90
index ea2b38645..f5e8d3795 100644
--- a/src/phase_damage_anisoductile.f90
+++ b/src/phase_damage_anisoductile.f90
@@ -4,7 +4,7 @@
 !> @brief material subroutine incorporating anisotropic ductile damage source mechanism
 !> @details to be done
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:constitutive_damage)  source_damage_anisoDuctile
+submodule(phase:damagee) anisoductile
 
   integer,                       dimension(:),           allocatable :: &
     source_damage_anisoDuctile_offset, &                                                            !< which source is my current damage mechanism?
@@ -184,4 +184,4 @@ module subroutine anisoductile_results(phase,group)
 
 end subroutine anisoductile_results
 
-end submodule source_damage_anisoDuctile
+end submodule anisoductile
diff --git a/src/phase_damage_isobrittle.f90 b/src/phase_damage_isobrittle.f90
index 0f192e2ae..5a01a425d 100644
--- a/src/phase_damage_isobrittle.f90
+++ b/src/phase_damage_isobrittle.f90
@@ -4,7 +4,7 @@
 !> @brief material subroutine incoprorating isotropic brittle damage source mechanism
 !> @details to be done
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:constitutive_damage) source_damage_isoBrittle
+submodule(phase:damagee) isobrittle
 
   integer,                       dimension(:),           allocatable :: &
     source_damage_isoBrittle_offset, &
@@ -180,4 +180,4 @@ module subroutine isobrittle_results(phase,group)
 
 end subroutine isobrittle_results
 
-end submodule source_damage_isoBrittle
+end submodule isobrittle
diff --git a/src/phase_damage_isoductile.f90 b/src/phase_damage_isoductile.f90
index 458b82e16..fe5ca93ed 100644
--- a/src/phase_damage_isoductile.f90
+++ b/src/phase_damage_isoductile.f90
@@ -4,7 +4,7 @@
 !> @brief material subroutine incorporating isotropic ductile damage source mechanism
 !> @details to be done
 !--------------------------------------------------------------------------------------------------
-submodule (constitutive:constitutive_damage) source_damage_isoDuctile
+submodule(phase:damagee) isoductile
 
   integer,                       dimension(:),           allocatable :: &
     source_damage_isoDuctile_offset, &                                                              !< which source is my current damage mechanism?
@@ -175,4 +175,4 @@ module subroutine isoductile_results(phase,group)
 
 end subroutine isoductile_results
 
-end submodule source_damage_isoDuctile
+end submodule isoductile
diff --git a/src/phase_mechanics.f90 b/src/phase_mechanics.f90
index 1d9248dd6..3bca63886 100644
--- a/src/phase_mechanics.f90
+++ b/src/phase_mechanics.f90
@@ -1,7 +1,7 @@
 !----------------------------------------------------------------------------------------------------
 !> @brief internal microstructure state for all plasticity constitutive models
 !----------------------------------------------------------------------------------------------------
-submodule(constitutive) mechanics
+submodule(phase) mechanics
 
   enum, bind(c); enumerator :: &
     ELASTICITY_UNDEFINED_ID, &
@@ -1615,7 +1615,7 @@ subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, &
       case (KINEMATICS_thermal_expansion_ID) kinematicsType
         me = material_phaseMemberAt(co,ip,el)
         ph = material_phaseAt(co,el)
-        call kinematics_thermal_expansion_LiAndItsTangent(my_Li, my_dLi_dS, ph,me)
+        call thermalexpansion_LiAndItsTangent(my_Li, my_dLi_dS, ph,me)
       case default kinematicsType
         my_Li = 0.0_pReal
         my_dLi_dS = 0.0_pReal
diff --git a/src/phase_mechanics_eigendeformation.f90 b/src/phase_mechanics_eigendeformation.f90
index ddae8ec06..9a86945dd 100644
--- a/src/phase_mechanics_eigendeformation.f90
+++ b/src/phase_mechanics_eigendeformation.f90
@@ -1,2 +1,2 @@
-submodule(constitutive:mechanics) eigendeformation
+submodule(phase:mechanics) eigendeformation
 end submodule eigendeformation
diff --git a/src/phase_mechanics_eigendeformation_cleavageopening.f90 b/src/phase_mechanics_eigendeformation_cleavageopening.f90
index b0ea2be87..185aee520 100644
--- a/src/phase_mechanics_eigendeformation_cleavageopening.f90
+++ b/src/phase_mechanics_eigendeformation_cleavageopening.f90
@@ -4,7 +4,7 @@
 !> @brief material subroutine incorporating kinematics resulting from opening of cleavage planes
 !> @details to be done
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:eigendeformation) cleavageopening
+submodule(phase:eigendeformation) cleavageopening
 
   integer, dimension(:), allocatable :: kinematics_cleavage_opening_instance
 
diff --git a/src/phase_mechanics_eigendeformation_slipplaneopening.f90 b/src/phase_mechanics_eigendeformation_slipplaneopening.f90
index 1cab5e731..8b8b1dd29 100644
--- a/src/phase_mechanics_eigendeformation_slipplaneopening.f90
+++ b/src/phase_mechanics_eigendeformation_slipplaneopening.f90
@@ -4,7 +4,7 @@
 !> @brief material subroutine incorporating kinematics resulting from opening of slip planes
 !> @details to be done
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:eigendeformation) slipplaneopening
+submodule(phase:eigendeformation) slipplaneopening
 
   integer, dimension(:), allocatable :: kinematics_slipplane_opening_instance
 
diff --git a/src/phase_mechanics_eigendeformation_thermalexpansion.f90 b/src/phase_mechanics_eigendeformation_thermalexpansion.f90
index 32bf30e75..5380330f2 100644
--- a/src/phase_mechanics_eigendeformation_thermalexpansion.f90
+++ b/src/phase_mechanics_eigendeformation_thermalexpansion.f90
@@ -3,7 +3,10 @@
 !> @brief material subroutine incorporating kinematics resulting from thermal expansion
 !> @details to be done
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:thermal) thermalexpansion
+submodule(phase:thermal) thermalexpansion
+  use prec
+  use YAML_types
+  use config
 
   integer, dimension(:), allocatable :: kinematics_thermal_expansion_instance
 
@@ -84,7 +87,7 @@ end function kinematics_thermal_expansion_init
 !--------------------------------------------------------------------------------------------------
 !> @brief constitutive equation for calculating the velocity gradient
 !--------------------------------------------------------------------------------------------------
-module subroutine kinematics_thermal_expansion_LiAndItsTangent(Li, dLi_dTstar, ph,me)
+module subroutine thermalexpansion_LiAndItsTangent(Li, dLi_dTstar, ph,me)
 
   integer, intent(in) :: ph, me
   real(pReal),   intent(out), dimension(3,3) :: &
@@ -92,9 +95,6 @@ module subroutine kinematics_thermal_expansion_LiAndItsTangent(Li, dLi_dTstar, p
   real(pReal),   intent(out), dimension(3,3,3,3) :: &
     dLi_dTstar                                                                                      !< derivative of Li with respect to Tstar (4th-order tensor defined to be zero)
 
-  integer :: &
-    phase, &
-    homog
   real(pReal) :: T, dot_T
 
   T     = current(ph)%T(me)
@@ -114,6 +114,6 @@ module subroutine kinematics_thermal_expansion_LiAndItsTangent(Li, dLi_dTstar, p
   end associate
   dLi_dTstar = 0.0_pReal
 
-end subroutine kinematics_thermal_expansion_LiAndItsTangent
+end subroutine thermalexpansion_LiAndItsTangent
 
 end submodule thermalexpansion
diff --git a/src/phase_mechanics_plastic.f90 b/src/phase_mechanics_plastic.f90
index 271c73391..035b03171 100644
--- a/src/phase_mechanics_plastic.f90
+++ b/src/phase_mechanics_plastic.f90
@@ -1,4 +1,4 @@
-submodule(constitutive:mechanics) plastic
+submodule(phase:mechanics) plastic
 
   interface
 
diff --git a/src/phase_mechanics_plastic_dislotungsten.f90 b/src/phase_mechanics_plastic_dislotungsten.f90
index 33379aae2..c57669edb 100644
--- a/src/phase_mechanics_plastic_dislotungsten.f90
+++ b/src/phase_mechanics_plastic_dislotungsten.f90
@@ -5,7 +5,7 @@
 !> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH
 !> @brief crystal plasticity model for bcc metals, especially Tungsten
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:plastic) dislotungsten
+submodule(phase:plastic) dislotungsten
 
   real(pReal), parameter :: &
     kB = 1.38e-23_pReal                                                                             !< Boltzmann constant in J/Kelvin
diff --git a/src/phase_mechanics_plastic_dislotwin.f90 b/src/phase_mechanics_plastic_dislotwin.f90
index db15e73b0..4889b7698 100644
--- a/src/phase_mechanics_plastic_dislotwin.f90
+++ b/src/phase_mechanics_plastic_dislotwin.f90
@@ -7,7 +7,7 @@
 !> @brief material subroutine incoprorating dislocation and twinning physics
 !> @details to be done
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:plastic) dislotwin
+submodule(phase:plastic) dislotwin
 
   real(pReal), parameter :: &
     kB = 1.38e-23_pReal                                                                             !< Boltzmann constant in J/Kelvin
diff --git a/src/phase_mechanics_plastic_isotropic.f90 b/src/phase_mechanics_plastic_isotropic.f90
index b79c14607..080a56cfd 100644
--- a/src/phase_mechanics_plastic_isotropic.f90
+++ b/src/phase_mechanics_plastic_isotropic.f90
@@ -7,7 +7,7 @@
 !! resolving the stress on the slip systems. Will give the response of phenopowerlaw for an
 !! untextured polycrystal
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:plastic) isotropic
+submodule(phase:plastic) isotropic
 
   type :: tParameters
     real(pReal) :: &
diff --git a/src/phase_mechanics_plastic_kinehardening.f90 b/src/phase_mechanics_plastic_kinehardening.f90
index ce6ea1935..de80a646d 100644
--- a/src/phase_mechanics_plastic_kinehardening.f90
+++ b/src/phase_mechanics_plastic_kinehardening.f90
@@ -5,7 +5,7 @@
 !> @brief  Phenomenological crystal plasticity using a power law formulation for the shear rates
 !! and a Voce-type kinematic hardening rule
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:plastic) kinehardening
+submodule(phase:plastic) kinehardening
 
   type :: tParameters
     real(pReal) :: &
diff --git a/src/phase_mechanics_plastic_none.f90 b/src/phase_mechanics_plastic_none.f90
index b09b4fc39..5c79c1b1e 100644
--- a/src/phase_mechanics_plastic_none.f90
+++ b/src/phase_mechanics_plastic_none.f90
@@ -4,7 +4,7 @@
 !> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH
 !> @brief Dummy plasticity for purely elastic material
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:plastic) none
+submodule(phase:plastic) none
 
 contains
 
diff --git a/src/phase_mechanics_plastic_nonlocal.f90 b/src/phase_mechanics_plastic_nonlocal.f90
index 4b44ee019..884bc6bc0 100644
--- a/src/phase_mechanics_plastic_nonlocal.f90
+++ b/src/phase_mechanics_plastic_nonlocal.f90
@@ -4,7 +4,7 @@
 !> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH
 !> @brief material subroutine for plasticity including dislocation flux
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:plastic) nonlocal
+submodule(phase:plastic) nonlocal
   use geometry_plastic_nonlocal, only: &
     nIPneighbors    => geometry_plastic_nonlocal_nIPneighbors, &
     IPneighborhood  => geometry_plastic_nonlocal_IPneighborhood, &
diff --git a/src/phase_mechanics_plastic_phenopowerlaw.f90 b/src/phase_mechanics_plastic_phenopowerlaw.f90
index a0318d266..a9538aeec 100644
--- a/src/phase_mechanics_plastic_phenopowerlaw.f90
+++ b/src/phase_mechanics_plastic_phenopowerlaw.f90
@@ -4,7 +4,7 @@
 !> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH
 !> @brief  phenomenological crystal plasticity formulation using a powerlaw fitting
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:plastic) phenopowerlaw
+submodule(phase:plastic) phenopowerlaw
 
   type :: tParameters
     real(pReal) :: &
diff --git a/src/phase_thermal.f90 b/src/phase_thermal.f90
index c0a999bde..5a7c7bc62 100644
--- a/src/phase_thermal.f90
+++ b/src/phase_thermal.f90
@@ -1,7 +1,7 @@
 !----------------------------------------------------------------------------------------------------
 !> @brief internal microstructure state for all thermal sources and kinematics constitutive models
 !----------------------------------------------------------------------------------------------------
-submodule(constitutive) thermal
+submodule(phase) thermal
 
   enum, bind(c); enumerator :: &
     THERMAL_UNDEFINED_ID ,&
diff --git a/src/phase_thermal_dissipation.f90 b/src/phase_thermal_dissipation.f90
index 3b4ddcfd1..ab2a13c64 100644
--- a/src/phase_thermal_dissipation.f90
+++ b/src/phase_thermal_dissipation.f90
@@ -4,7 +4,7 @@
 !> @brief material subroutine for thermal source due to plastic dissipation
 !> @details to be done
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:thermal) dissipation
+submodule(phase:thermal) dissipation
 
   integer,           dimension(:),   allocatable :: &
     source_thermal_dissipation_offset, &                                                            !< which source is my current thermal dissipation mechanism?
diff --git a/src/phase_thermal_externalheat.f90 b/src/phase_thermal_externalheat.f90
index 054b91e01..5707efb11 100644
--- a/src/phase_thermal_externalheat.f90
+++ b/src/phase_thermal_externalheat.f90
@@ -4,7 +4,7 @@
 !> @author Philip Eisenlohr, Michigan State University
 !> @brief material subroutine for variable heat source
 !--------------------------------------------------------------------------------------------------
-submodule(constitutive:thermal) externalheat
+submodule(phase:thermal) externalheat
 
 
   integer,           dimension(:),   allocatable :: &