Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 4 additions & 4 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -41,11 +41,11 @@ produces output that includes actual program syntax:
```fortran
f = (x**2)/2 ! <-- scalar function
v = x ! <-- vector function
.SSS. (v .dot. .grad. f) * dV = .3333333330205934
.SSS. ( f * .div. v) * dV = .16666666739857125
-.SS. (f .x. (v .dot. dA)) = -.5000000004191649
.SSS. (v .dot. .grad. f) * dV = .33333322
.SSS. ( f * .div. v) * dV = .16666671
-.SS. (f .x. (v .dot. dA)) = -.4999999
----------------------------------------------------
sum = -.2220446049250313E-15 (residual)
sum = .29802323E-7 (residual)
```
where the small residual evidences a highly accurate approximation.

Expand Down
22 changes: 11 additions & 11 deletions doc/uml/api.md
Original file line number Diff line number Diff line change
Expand Up @@ -12,16 +12,16 @@ vector_3D_t <|-- gradient_3D_t : is a
class scalar_3D_t{
operator(.grad.) (scalar_3D_t) gradient_3D_t
operator(+) (scalar_3D_t, scalar_3D_t) scalar_3D_t
operator(*) (double precision, scalar_3D_t) scalar_3D_t
operator(*) (scalar_3D_t, double precision) scalar_3D_t
operator(*) (real, scalar_3D_t) scalar_3D_t
operator(*) (scalar_3D_t, real) scalar_3D_t
operator(*) (integer, scalar_3D_t) scalar_3D_t
operator(*) (scalar_3D_t, integer) scalar_3D_t
values() double precision(:,:,:)
grid(component : integer, coordinate : integer) double precision(:)
values() real(:,:,:)
grid(component : integer, coordinate : integer) real(:)
consistent() : logical
to_faces(direction : integer) double precision(:,:,:)
to_faces(direction : integer) real(:,:,:)
to_file(name : character(len=:)) file_t
scalar_3D_t(initializer : scalar_3D_initializer_i, order : integer, cells : integer(:), x_min : double precision(:), x_max : double precision(:)) scalar_3D_t
scalar_3D_t(initializer : scalar_3D_initializer_i, order : integer, cells : integer(:), x_min : real(:), x_max : real(:)) scalar_3D_t
scalar_3D_t(initializer : scalar_3D_initializer_i, mold : scalar_3D_t) scalar_3D_t
}

Expand All @@ -30,18 +30,18 @@ class vector_3D_t{
operator(.dot.) (vector_3D_t, vector_3D_t) : scalar_3D_t
operator(*)(scalar_3D_t, vector_3D_t) vector_3D_t
operator(*)(vector_3D_t, scalar_3D_t) vector_3D_t
grid(component : integer, coordinate : integer) double precision(:)
values(direction : integer) double precision(:,:,:)
grid(component : integer, coordinate : integer) real(:)
values(direction : integer) real(:,:,:)
consistent() : logical
to_centers_extended() double precision(:,:,:,:)
to_centers_extended() real(:,:,:,:)
to_file(name : character(len=:)) file_t
vector_3D_t(initializer : vector_3D_initializer_i, order : integer, cells : integer(:), x_min : double precision(:), x_max : double precision(:))
vector_3D_t(initializer : vector_3D_initializer_i, order : integer, cells : integer(:), x_min : real(:), x_max : real(:))
vector_3D_t(initializer : vector_3D_initializer_i, mold : vector_3D_t)
vector_3D_t(initializer : vector_3D_initializer_i, mold : scalar_3D_t)
}

class divergence_3D_t{
assignment(=)(scalar_3D_t, divergence_3D_t)
divergence_3D_t(initializer : vector_3D_initializer_i, order : integer, cells : integer(:), x_min : double precision(:), x_max : double precision(:))
divergence_3D_t(initializer : vector_3D_initializer_i, order : integer, cells : integer(:), x_min : real(:), x_max : real(:))
divergence_3D_t(initializer : vector_3D_initializer_i, mold : vector_3D_t)
}
16 changes: 8 additions & 8 deletions doc/uml/class-diagram.md
Original file line number Diff line number Diff line change
Expand Up @@ -45,21 +45,21 @@ class vector_1D_t{
}

class gradient_1D_t{
- weights : double precision[]
- weights : real[]
}

class mimetic_matrix_1D_t{
- upper_ :: double precision[]
- inner_ :: double precision[]
- lower_ :: double precision[]
- upper_ :: real[]
- inner_ :: real[]
- lower_ :: real[]
}

class gradient_operator_1D_t{
+ operator(.x.) double precision[]
+ assemble() double precision[] "2D array"
+ operator(.x.) real[]
+ assemble() real[] "2D array"
}

class divergence_operator_1D_t{
+ operator(.x.) double precision[]
+ assemble() double precision[] "2D array"
+ operator(.x.) real[]
+ assemble() real[] "2D array"
}
18 changes: 9 additions & 9 deletions example/2D-advection-diffusion.F90
Original file line number Diff line number Diff line change
Expand Up @@ -9,10 +9,10 @@ module fields_m
contains

pure function scalar_field(x,y) result(gaussian)
double precision, intent(in) :: x(:), y(:)
double precision gaussian(size(x),size(y))
double precision, parameter :: pi = acos(-1D0)
double precision, parameter :: x0 = -pi/2, y0 = -pi/2, sigma = pi/8
real, intent(in) :: x(:), y(:)
real gaussian(size(x),size(y))
real, parameter :: pi = acos(-1E0)
real, parameter :: x0 = -pi/2, y0 = -pi/2, sigma = pi/8
do concurrent(integer :: j=1:size(y)) default(none) shared(x,y,gaussian)
associate(r => sqrt((x-x0)**2 + (y(j)-y0)**2))
gaussian(:,j) = exp(-(r**2)/(2*sigma**2))
Expand All @@ -21,8 +21,8 @@ pure function scalar_field(x,y) result(gaussian)
end function

pure function taylor_green_velocity(x,y) result(velocity)
double precision, intent(in) :: x(:), y(:)
double precision velocity(size(x),size(y),space_dimension)
real, intent(in) :: x(:), y(:)
real velocity(size(x),size(y),space_dimension)
do concurrent(integer :: i=1:size(x), j=1:size(y))
velocity(i,j,:) = [10*sin(x(i))*cos(y(j)), -10*cos(x(i))*sin(y(j))]
end do
Expand All @@ -41,7 +41,7 @@ program advection_diffusion_2D
procedure(scalar_2D_initializer_i), pointer :: scalar_2D_initializer
procedure(vector_2D_initializer_i), pointer :: velocity_2D_initializer
type(scalar_2D_t) s
double precision, parameter :: pi = acos(-1D0)
real, parameter :: pi = acos(-1E0)

scalar_2D_initializer => scalar_field
velocity_2D_initializer => taylor_green_velocity
Expand All @@ -60,7 +60,7 @@ program advection_diffusion_2D

advance_time: &
block
double precision :: dt = 1D-4
real :: dt = 1E-4
integer step

do step = 1, 500
Expand Down Expand Up @@ -88,7 +88,7 @@ pure function d_dt(s, v) result(ds_dt)
type(scalar_2D_t), intent(in) :: s
type(vector_2D_t), intent(in) :: v
type(scalar_2D_t) ds_dt
double precision, parameter :: D = 0.5D0
real, parameter :: D = 0.5E0
ds_dt = .div. (D * .grad. s) - .div. (v * s)
end function

Expand Down
20 changes: 10 additions & 10 deletions example/2D-sink.F90
Original file line number Diff line number Diff line change
Expand Up @@ -8,27 +8,27 @@ module sink_2D_functions_m
implicit none

integer, parameter :: space_dimension = 2
double precision, parameter :: pi = acos(-1D0)
real, parameter :: pi = acos(-1E0)

contains

pure function velocity(x,y) result(v)
double precision, intent(in) :: x(:), y(:)
double precision v(size(x),size(y),space_dimension)
double precision, parameter :: Q = 1D0
real, intent(in) :: x(:), y(:)
real v(size(x),size(y),space_dimension)
real, parameter :: Q = 1E0
do concurrent(integer :: i=1:size(x), j=1:size(y))
associate(r => sqrt(x(i)**2 + y(j)**2))
call_julienne_assert(r /= 0D0)
call_julienne_assert(r /= 0E0)
v(i,j,:) = -(Q/(2*pi))*[x(i), y(j)]/(x(i)**2 + y(j)**2)
end associate
end do
end function

pure function divergence(x,y) result(div_v)
double precision, intent(in) :: x(:), y(:)
double precision div_v(size(x),size(y))
call_julienne_assert(.not. any(x == 0D0 .and. y == 0D0))
div_v = 0D0
real, intent(in) :: x(:), y(:)
real div_v(size(x),size(y))
call_julienne_assert(.not. any(x == 0E0 .and. y == 0E0))
div_v = 0E0
end function

end module
Expand All @@ -45,7 +45,7 @@ program sink_2D
divergence_2D_initializer => divergence
vector_2D_initializer => velocity

associate(v => vector_2D_t(vector_2D_initializer, order=order, cells=[11,11], x_min=[-1D0,-1D0], x_max=[1D0,1D0]))
associate(v => vector_2D_t(vector_2D_initializer, order=order, cells=[11,11], x_min=[-1E0,-1E0], x_max=[1E0,1E0]))
associate(div_v => .div. v, expected_divergence => divergence_2D_t(divergence_2D_initializer, mold=v))
associate(v_file => v%to_file("v"),div_v_file => div_v%to_file(".div. v"), expected_divergence_file => expected_divergence%to_file("expected .div. v"))
call v_file%write_lines("example/scripts/sink-velocity.csv")
Expand Down
10 changes: 5 additions & 5 deletions example/2D-stagnation-point.F90
Original file line number Diff line number Diff line change
Expand Up @@ -9,16 +9,16 @@ module velocity_potential_m
contains

pure function potential(x,y) result(phi)
double precision, intent(in) :: x(:), y(:)
double precision phi(size(x),size(y))
real, intent(in) :: x(:), y(:)
real phi(size(x),size(y))
do concurrent(integer :: j=1:size(y)) default(none) shared(x,y,phi)
phi(:,j) = (x**2 - y(j)**2)/2
end do
end function

pure function potential_gradient(x,y) result(grad_phi)
double precision, intent(in) :: x(:), y(:)
double precision grad_phi(size(x),size(y),space_dimension)
real, intent(in) :: x(:), y(:)
real grad_phi(size(x),size(y),space_dimension)
do concurrent(integer :: i=1:size(x), j=1:size(y))
grad_phi(i,j,:) = [x(i), -y(j)]
end do
Expand All @@ -38,7 +38,7 @@ program stagnation_point_2D
vector_2D_initializer => potential_gradient

associate( &
phi => scalar_2D_t(scalar_2D_initializer, order=4, cells=[20,20], x_min=[-2D0,-2D0], x_max=[2D0,2D0]) &
phi => scalar_2D_t(scalar_2D_initializer, order=4, cells=[20,20], x_min=[-2E0,-2E0], x_max=[2E0,2E0]) &
)
associate( &
v => .grad. phi &
Expand Down
22 changes: 11 additions & 11 deletions example/burgers-1D.F90
Original file line number Diff line number Diff line change
Expand Up @@ -10,8 +10,8 @@ module initial_condition_m

pure function initial_condition(x)
!! Initial solution to Burgers equation
double precision, intent(in) :: x(:)
double precision, allocatable :: initial_condition(:)
real, intent(in) :: x(:)
real, allocatable :: initial_condition(:)
initial_condition = 10*sin(x)
! To change this function, please edit only the right-hand-side (RHS) expression,
! keeping the rest in place for proper display of the function at runtime.
Expand Down Expand Up @@ -66,14 +66,14 @@ program burgers_1D

block
procedure(scalar_1D_initializer_i), pointer :: scalar_1D_initializer
double precision, parameter :: pi = acos(-1D0), nu=1D0, t_final=0.6D0
double precision, allocatable :: u_surface(:,:)
double precision dt
real, parameter :: pi = acos(-1E0), nu=1E0, t_final=0.6E0
real, allocatable :: u_surface(:,:)
real dt
type(scalar_1D_t) u
integer step, n

scalar_1D_initializer => initial_condition
u = scalar_1D_t(scalar_1D_initializer, order, x_min=0D0, x_max=2*pi, cells=199)
u = scalar_1D_t(scalar_1D_initializer, order, x_min=0E0, x_max=2*pi, cells=199)
dt = diffusion_stability_limit(nu, u%dx(), order)

associate(steps => ceiling(t_final/dt))
Expand Down Expand Up @@ -148,17 +148,17 @@ program burgers_1D

pure function d_dt(u, nu) result(du_dt)
type(scalar_1D_t), intent(in) :: u
double precision, intent(in) :: nu
real, intent(in) :: nu
type(scalar_1D_t) du_dt
du_dt = nu*d2_dx2(u) - d_dx((u**2)/2)
end function

pure function diffusion_stability_limit(diffusivity,delta_x,order_of_accuracy) result(stable_time_step)
double precision, intent(in) :: diffusivity, delta_x
real, intent(in) :: diffusivity, delta_x
integer, intent(in) :: order_of_accuracy
double precision stable_time_step
double precision, parameter, dimension(*) :: stability_limit=[2.,2.,2.5,2.79] ! third value needs to be checked
double precision, parameter :: safety_factor = 0.9
real stable_time_step
real, parameter, dimension(*) :: stability_limit=[2.,2.,2.5,2.79] ! third value needs to be checked
real, parameter :: safety_factor = 0.9
! See Moin, P. (2010) Fundamentals of Engineering Numerical Analysis, 2nd ed., pp. 111-116.
stable_time_step = safety_factor*stability_limit(order_of_accuracy)*(delta_x**2)/(4*diffusivity)
end function
Expand Down
16 changes: 8 additions & 8 deletions example/div-grad-laplacian-1D.F90
Original file line number Diff line number Diff line change
Expand Up @@ -8,25 +8,25 @@ module functions_m

pure function f(x)
!! Function for differentiation
double precision, intent(in) :: x(:)
double precision, allocatable :: f(:)
real, intent(in) :: x(:)
real, allocatable :: f(:)
f = (x**3)/6 + (x**2)/2 + 1
! To change functions, edit only the right-hand-side (RHS) expression.
! Please keep the rest place for proper display of the function at runtime.
end function

double precision elemental function df_dx(x)
real elemental function df_dx(x)
!! 1st-derivative function
double precision, intent(in) :: x
real, intent(in) :: x
df_dx = (x**2)/2 + x
! To change derivative functions, edit only the RHS.
! Please keep the rest place for proper display of the function at runtime.
! Also, ensure the new RHS expresses the first derivative of f above.
end function

double precision elemental function d2f_dx2(x)
real elemental function d2f_dx2(x)
!! 2nd-derivative function
double precision, intent(in) :: x
real, intent(in) :: x
d2f_dx2 = x + 1
! To edit the above function, edit only the right-hand side expression,
! Please keep the rest place for proper display of the function at runtime.
Expand Down Expand Up @@ -93,7 +93,7 @@ program div_grad_laplacian_1D
subroutine output(order)
integer, intent(in) :: order

associate( s => scalar_1D_t(scalar_1D_initializer, order=order, cells=20, x_min=0D0, x_max=20D0))
associate( s => scalar_1D_t(scalar_1D_initializer, order=order, cells=20, x_min=0E0, x_max=20E0))
associate( grad_s => .grad. s &
,laplacian_s => .laplacian. s)
associate( s_grid => s%grid() &
Expand Down Expand Up @@ -136,7 +136,7 @@ subroutine output(order)
end subroutine

pure function tabulate(headings, abscissa, expected, actual) result(file)
double precision, intent(in), dimension(:) :: abscissa, expected, actual
real, intent(in), dimension(:) :: abscissa, expected, actual
type(string_t), intent(in) :: headings(:)
type(file_t) file
integer line
Expand Down
14 changes: 7 additions & 7 deletions example/extended-gauss-divergence.F90
Original file line number Diff line number Diff line change
Expand Up @@ -6,14 +6,14 @@ module integrand_operands_m
contains

pure function scalar(x) result(f)
double precision, intent(in) :: x(:)
double precision, allocatable :: f(:)
real, intent(in) :: x(:)
real, allocatable :: f(:)
f = (x**2)/2 ! <-- scalar function
end function

pure function vector(x) result(v)
double precision, intent(in) :: x(:)
double precision, allocatable :: v(:)
real, intent(in) :: x(:)
real, allocatable :: v(:)
v = x ! <-- vector function
end function

Expand All @@ -36,15 +36,15 @@ program extended_gauss_divergence
!! Define default initializations that can be overridden with the command-line arguments
!! detailed by the usage information below
integer :: cells_=200, order_=4
double precision :: x_min_=0D0, x_max_=1D0
real :: x_min_=0E0, x_max_=1E0
end type

type text_flags_t
logical div_, grad_, vf_
end type

type(command_line_t) command_line
double precision SSS_v_dot_grad_f_dV, SSS_f_div_v_dV, SS_f_v_dot_dA
real SSS_v_dot_grad_f_dV, SSS_f_div_v_dV, SS_f_v_dot_dA

if (command_line%argument_present([character(len=len("--help")) :: ("--help"), "-h"])) then
stop new_line('') // new_line('') &
Expand All @@ -53,7 +53,7 @@ program extended_gauss_divergence
// ' --example extended-gauss-divergence \' // new_line('') &
// ' --compiler flang-new \' // new_line('') &
// ' --flag "-O3" \' // new_line('') &
// ' -- [--help|-h] | [[--cells <integer>] [--order <integer>] [--xmin <double precision>] [--xmax <double precision>] [--div|d] [--grad|g] [--vf|f]]' &
// ' -- [--help|-h] | [[--cells <integer>] [--order <integer>] [--xmin <real>] [--xmax <real>] [--div|d] [--grad|g] [--vf|f]]' &
// new_line('') // new_line('') &
// 'where pipes (|) separate square-bracketed optional arguments and angular brackets indicate user input values.' // new_line('')
end if
Expand Down
Loading