Data type refactor

configure.ac

@@ -474,6 +474,8 @@ AC_CONFIG_FILES([tests/unit/stack/Makefile\
         tests/unit/octree/Makefile\
         tests/unit/vector/Makefile\
         tests/unit/matrix/Makefile\
+        tests/unit/tensor3/Makefile\
+        tests/unit/tensor4/Makefile\
         tests/unit/case_file_utils/Makefile\
         tests/unit/scratch_registry/Makefile\
         tests/unit/registry/Makefile\

contrib/calc_lift_from_field/calc_lift_from_field.f90

@@ -142,7 +142,6 @@ program calc_lift_from_field
      call copy(v%x,field_data%v%x,dof%size())
      call copy(w%x,field_data%w%x,dof%size())
      call copy(p%x,field_data%p%x,dof%size())
-     drag_torq = 0.0_rp
      !set coords to somoething big
      drag_torq%x(:,2) = 1e15
 

src/Makefile.am

@@ -36,8 +36,13 @@ neko_fortran_SOURCES = \
 	common/system.f90\
 	common/runtime_statistics.f90\
 	common/user_access_singleton.f90\
+	data_types/array.f90\
 	data_types/host_array.f90\
 	data_types/device_array.f90\
+	data_types/vector.f90\
+	data_types/matrix.f90\
+	data_types/tensor3.f90\
+	data_types/tensor4.f90\
 	math/mxm_wrapper.F90\
 	sem/speclib.f90 \
 	qoi/drag_torque.f90\
@@ -117,6 +122,7 @@ neko_fortran_SOURCES = \
 	math/bcknd/cpu/cpu_opgrad.f90\
 	math/bcknd/cpu/convect_scalar.f90\
 	math/bcknd/cpu/set_convect_rst.f90\
+	math/bcknd/cpu/cpu_matrix_math.f90\
 	math/bcknd/sx/sx_cdtp.f90\
 	math/bcknd/sx/sx_conv1.f90\
 	math/bcknd/sx/sx_dudxyz.f90\
@@ -139,15 +145,13 @@ neko_fortran_SOURCES = \
 	math/bcknd/sx/fdm_sx.f90\
 	math/bcknd/xsmm/fdm_xsmm.f90\
 	math/schwarz.f90\
-	math/vector.f90\
 	common/particles.f90\
 	math/vector_series.f90\
 	math/vector_list.f90\
 	math/vector_math.f90\
 	registries/registry.f90\
 	registries/registry_entry.f90\
 	registries/scratch_registry.f90\
-	math/matrix.f90\
 	math/matrix_math.f90\
 	math/signed_distance.f90\
 	multigrid/tree_amg.f90\

src/bc/neumann.f90

@@ -43,6 +43,7 @@ module neumann
   use vector, only : vector_t
   use neko_config, only : NEKO_BCKND_DEVICE
   use device_math, only : device_cfill, device_copy
+  use vector_math, only : vector_cfill, vector_copy
   use device, only : device_memcpy, DEVICE_TO_HOST
   use device_neumann, only : device_neumann_apply_scalar, &
        device_neumann_apply_vector
@@ -360,7 +361,7 @@ subroutine neumann_finalize(this)
     ! Allocate flux vectors and assign to initial constant values
     do i = 1, size(this%init_flux_)
        call this%flux(i)%init(this%facet_node_msk(0))
-       this%flux(i) = this%init_flux_(i)
+       call vector_cfill(this%flux(i), this%init_flux_(i))
     end do
 
     this%uniform_0 = .true.
@@ -383,7 +384,7 @@ subroutine neumann_set_flux_scalar(this, flux, comp)
             "neumann_set_flux_scalar")
     end if
 
-    this%flux(comp) = flux
+    call vector_cfill(this%flux(comp), flux)
     ! If we were uniform zero before, and this comp is set to zero, we are still
     ! uniform zero
     this%uniform_0 = abscmp(flux, 0.0_rp) .and. this%uniform_0
@@ -404,7 +405,7 @@ subroutine neumann_set_flux_array(this, flux, comp)
             "neuman_set_flux_array")
     end if
 
-    this%flux(comp) = flux
+    call vector_copy(this%flux(comp), flux)
 
     ! Once a flux is set explicitly, we no longer assume it is uniform zero.
     this%uniform_0 = .false.

src/bc/non_normal_aligned.f90

@@ -44,7 +44,7 @@ module non_normal_aligned
        json_get_or_lookup
   use utils, only : neko_error
   use vector, only : vector_t
-  use vector_math, only : vector_masked_gather_copy_0
+  use vector_math, only : vector_cfill, vector_masked_gather_copy_0
   use scratch_registry, only : neko_scratch_registry
   use field, only : field_t
   use global_interpolation, only : GLOB_INTERP_TOL, GLOB_INTERP_PAD
@@ -421,9 +421,9 @@ subroutine non_normal_aligned_finalize(this)
     call this%value_z%init(this%bc_z%msk(0), 'non_normal_aligned_z')
 
     if (this%use_constant_value) then
-       this%value_x = this%constant_value(1)
-       this%value_y = this%constant_value(2)
-       this%value_z = this%constant_value(3)
+       call vector_cfill(this%value_x, this%constant_value(1))
+       call vector_cfill(this%value_y, this%constant_value(2))
+       call vector_cfill(this%value_z, this%constant_value(3))
     else if (this%read_values_from_field) then
        call neko_scratch_registry%request_field(value_x_field, scratch_idx(1), &
             .true.)

src/bc/vector_bc_projector.f90

@@ -980,10 +980,6 @@ subroutine coupled_vector_bc_projector_rebuild_basis(this)
     call rzero(normal_y_field%x, dof_size)
     call rzero(normal_z_field%x, dof_size)
 
-    this%n = 0.0_rp
-    this%t1 = 0.0_rp
-    this%t2 = 0.0_rp
-
     do i = 1, this%bcs%size()
        bc => this%bcs%get(i)
 

src/common/particles.f90

@@ -35,6 +35,7 @@ module particles
   use num_types, only : rp
   use vector, only : vector_t
   use utils, only : neko_error
+  use device, only : HOST_TO_DEVICE
   implicit none
   private
 
@@ -112,31 +113,34 @@ subroutine particles_init(this, x, y, z, time_order, u, v, w, diameter, &
     call this%acc_zlaglag%init(this%n)
     call this%d%init(this%n)
     call this%rho%init(this%n)
-    this%x = x
-    this%y = y
-    this%z = z
+
+    this%x%x = x
+    this%y%x = y
+    this%z%x = z
     if (present(u) .and. present(v) .and. present(w)) then
-       this%u = u
-       this%v = v
-       this%w = w
+       this%u%x = u
+       this%v%x = v
+       this%w%x = w
     else
-       this%u = 0.0_rp
-       this%v = 0.0_rp
-       this%w = 0.0_rp
+       this%u%x = 0.0_rp
+       this%v%x = 0.0_rp
+       this%w%x = 0.0_rp
     end if
     if (present(diameter)) then
-       this%d = diameter
+       this%d%x = diameter
     else
-       this%d = 0.0_rp
+       this%d%x = 0.0_rp
     end if
     if (present(density)) then
-       this%rho = density
+       this%rho%x = density
     else
-       this%rho = 0.0_rp
+       this%rho%x = 0.0_rp
     end if
     do i = 1, this%n
        this%ids(i) = i
     end do
+
+    call this%device_sync(HOST_TO_DEVICE)
   end subroutine particles_init
 
   subroutine particles_free(this)