Skip to content
Snippets Groups Projects
Verified Commit e22b0177 authored by Dorian Stoll's avatar Dorian Stoll
Browse files

tug: julia: Remove not fully implemented approaches

parent 24e5910b
No related branches found
No related tags found
No related merge requests found
...@@ -275,55 +275,6 @@ function ThomasAlgorithm( ...@@ -275,55 +275,6 @@ function ThomasAlgorithm(
return x_vec return x_vec
end end
# BTCS solution for 1D grid
function BTCS_1D!(grid::Grid{T}, bc::Boundary{T}, timestep::T, solver::SOLVER) where {T}
length::Int = grid_getLength(grid)
sx::T = timestep / (getGridDelta(grid) * getGridDelta(grid))
concentrations_t1 = Vector{T}(undef, length)
v = Vector{T}(undef, length)
alpha = getGridAlpha(grid)
bcLeft = getBoundarySide(bc, BC_SIDE_LEFT)
bcRight = getBoundarySide(bc, BC_SIDE_RIGHT)
concentrations::AbstractMatrix{T} = getGridConcentrations(grid)
rowIndex::Int = 1
solverFunc::Function
if solver == JULIA_BUILTIN_SOLVER
solverFunc = JuliaBuiltinAlgorithm
elseif solver == THOMAS_ALGORITHM_SOLVER
solverFunc = ThomasAlgorithm
else
error("Invalid solver type!")
end
A::AbstractMatrix{T} = createCoeffMatrix(alpha, bcLeft, bcRight, length, rowIndex, sx) # this is exactly same as in 2D
for i = 1:length
b[i] = concentrations[i, 0]
end
if getBoundaryElementType(bc, BC_SIDE_LEFT, 1) == BC_TYPE_CONSTANT
b[1] += 2 * sx * alpha[1, 1] * getBoundaryElementValue(bcLeft[1])
end
if getBoundaryElementType(bc, BC_SIDE_RIGHT, 1) == BC_TYPE_CONSTANT
b[length] += 2 * sx * alpha[length, 1] * getBoundaryElementValue(bcRight[1])
end
concentrations = solverFunc(A, b)
# for j = 1:length
# concentrations[j, 1] = concentrations_t1[j]
# end
# grid_setConcentrations!(grid, concentrations)
end
# BTCS solution for 2D grid # BTCS solution for 2D grid
function BTCS_2D!(grid::Grid{T}, bc::Boundary{T}, timestep::T, solver::SOLVER) where {T} function BTCS_2D!(grid::Grid{T}, bc::Boundary{T}, timestep::T, solver::SOLVER) where {T}
rowMax::Int = getGridRow(grid) rowMax::Int = getGridRow(grid)
...@@ -412,20 +363,9 @@ end ...@@ -412,20 +363,9 @@ end
# differentiate between 1D and 2D grid # differentiate between 1D and 2D grid
function BTCS!(grid::Grid{T}, bc::Boundary{T}, timestep::T, solver::SOLVER) where {T} function BTCS!(grid::Grid{T}, bc::Boundary{T}, timestep::T, solver::SOLVER) where {T}
if getGridDim(grid) == 1 if getGridDim(grid) == 2
BTCS_1D!(grid, bc, timestep, solver)
elseif getGridDim(grid) == 2
BTCS_2D!(grid, bc, timestep, solver) BTCS_2D!(grid, bc, timestep, solver)
else else
error("Error: Only 1- and 2-dimensional grids are defined!") error("Error: Only 2-dimensional grids are defined!")
end end
end end
# entry point for Julia solver
function BTCS_Julia!(grid::Grid{T}, bc::Boundary{T}, timestep::T) where {T}
BTCS!(grid, bc, timestep, JULIA_BUILTIN_SOLVER)
end
function BTCS_Thomas!(grid::Grid{T}, bc::Boundary{T}, timestep::T) where {T}
BTCS!(grid, bc, timestep, THOMAS_ALGORITHM_SOLVER)
end
using Printf using Printf
# Enum defining the implemented solution approaches.
@enum APPROACH begin
FTCS_APPROACH # Forward Time-Centered Space
BTCS_APPROACH # Backward Time-Centered Space
CRANK_NICOLSON_APPROACH # Crank-Nicolson method
end
# Enum holding different options for .csv output. # Enum holding different options for .csv output.
@enum CSV_OUTPUT begin @enum CSV_OUTPUT begin
CSV_OUTPUT_OFF # do not produce csv output CSV_OUTPUT_OFF # do not produce csv output
...@@ -33,7 +26,6 @@ end ...@@ -33,7 +26,6 @@ end
# #
# T is the type of the internal data structures for grid, boundary condition and timestep # T is the type of the internal data structures for grid, boundary condition and timestep
mutable struct Simulation{T<:Real} mutable struct Simulation{T<:Real}
approach::APPROACH
solver::SOLVER solver::SOLVER
timestep::T timestep::T
...@@ -46,8 +38,6 @@ mutable struct Simulation{T<:Real} ...@@ -46,8 +38,6 @@ mutable struct Simulation{T<:Real}
grid::Grid{T} grid::Grid{T}
bc::Boundary{T} bc::Boundary{T}
approach_names::Dict{APPROACH,String}
# Set up a simulation environment. The timestep and number of iterations must be set. For the # Set up a simulation environment. The timestep and number of iterations must be set. For the
# BTCS approach, the Thomas algorithm is used as the default linear equation solver as this is # BTCS approach, the Thomas algorithm is used as the default linear equation solver as this is
# faster for tridiagonal matrices. CSV output, console output and time measure are off by default. # faster for tridiagonal matrices. CSV output, console output and time measure are off by default.
...@@ -55,17 +45,14 @@ mutable struct Simulation{T<:Real} ...@@ -55,17 +45,14 @@ mutable struct Simulation{T<:Real}
# Arguments: # Arguments:
# grid: Valid grid object # grid: Valid grid object
# bc: Valid boundary condition object # bc: Valid boundary condition object
# approach: Set the SLE scheme to be used
# solver: Set the solver to be used # solver: Set the solver to be used
function Simulation{T}( function Simulation{T}(
grid::Grid{T}, grid::Grid{T},
bc::Boundary{T}, bc::Boundary{T},
; ;
approach::APPROACH = BTCS_APPROACH,
solver::SOLVER = THOMAS_ALGORITHM_SOLVER, solver::SOLVER = THOMAS_ALGORITHM_SOLVER,
) where {T<:Real} ) where {T<:Real}
new( new(
approach,
solver, solver,
-1, -1,
-1, -1,
...@@ -75,11 +62,6 @@ mutable struct Simulation{T<:Real} ...@@ -75,11 +62,6 @@ mutable struct Simulation{T<:Real}
TIME_MEASURE_OFF, TIME_MEASURE_OFF,
grid, grid,
bc, bc,
Dict(
FTCS_APPROACH => "FTCS",
BTCS_APPROACH => "BTCS",
CRANK_NICOLSON_APPROACH => "CRNI",
),
) )
end end
end end
...@@ -115,65 +97,7 @@ function setSimulationTimestep!(this::Simulation{T}, timestep::T) where {T} ...@@ -115,65 +97,7 @@ function setSimulationTimestep!(this::Simulation{T}, timestep::T) where {T}
error("Timestep has to be greater than zero.") error("Timestep has to be greater than zero.")
end end
if this.approach == FTCS_APPROACH || this.approach == CRANK_NICOLSON_APPROACH this.timestep = timestep
cfl::T
maxAlpha::T
if getGridDim(this.grid) == 1
deltaSquare::T = getGridDelta(grid)
maxAlpha = findmax(getGridAlpha(grid))[1]
# Courant-Friedrichs-Lewy condition
cfl = deltaSquare / (4 * maxAlpha)
elseif grid_getDim(grid) == 2
deltaColSquare::T = getGridDeltaCol(grid) * getGridDeltaCol(grid)
# will be 0 if 1D, else ...
deltaRowSquare::T = getGridDeltaRow(grid) * getGridDeltaRow(grid)
minDeltaSquare::T = min(deltaColSquare, deltaRowSquare)
maxAlpha = max(findmax(getGridAlphaX(grid))[1], findmax(getGridAlphaY(grid))[1])
cfl = minDeltaSquare / (4 * maxAlpha)
end
dim::String = @sprintf "%dD" getGridDim(grid)
apporachPrefix::String = this.approach_names[this.approach]
println(stdout, apporachPrefix, "_", dim, " :: CFL condition: ", cfl)
println(stdout, approachPrefix, "_", dim, " :: required dt=", timestep)
if timestep > cfl
this.innerIterations = Int(ceil(timestep / cfl))
this.timestep = timestep / this.innerIterations
println(
stderr,
"Warning :: Timestep was adjusted, because of stability conditions. Time duration was approximately preserved by adjusting internal number of iterations.",
)
println(
stdout,
apporachPrefix,
"_",
dim,
" :: Required ",
this.innerIterations,
" inner iterations with dt=",
this.timestep,
)
else
this.timestep = timestep
println(
stdout,
apporachPrefix,
"_",
dim,
" :: No inner iterations required, dt=",
timestep,
)
end
else
this.timestep = timestep
end
end end
# Currently set time step is returned. # Currently set time step is returned.
...@@ -215,19 +139,16 @@ function createCSVFile(this::Simulation{T})::String where {T} ...@@ -215,19 +139,16 @@ function createCSVFile(this::Simulation{T})::String where {T}
appendIdent::Int = 0 appendIdent::Int = 0
appendIdentString::String = "" appendIdentString::String = ""
approachString::String = this.approach_names[this.approach]
row::String = string(getGridRow(this.grid)) row::String = string(getGridRow(this.grid))
col::String = string(getGridCol(this.grid)) col::String = string(getGridCol(this.grid))
numIterations::String = string(this.iterations) numIterations::String = string(this.iterations)
filename::String = approachString * "_" * row * "_" * col * "_" * numIterations * ".csv" filename::String = row * "_" * col * "_" * numIterations * ".csv"
while isfile(filename) while isfile(filename)
appendIdent += 1 appendIdent += 1
appendIdentString = string(appendIdent) appendIdentString = string(appendIdent)
filename = filename =
approachString *
"_" *
row * row *
"_" * "_" *
col * col *
...@@ -314,55 +235,16 @@ function runSimulation!(this::Simulation{T}) where {T} ...@@ -314,55 +235,16 @@ function runSimulation!(this::Simulation{T}) where {T}
filename = createCSVFile(this) filename = createCSVFile(this)
end end
if this.approach == FTCS_APPROACH for i = 1:this.iterations
for i = 1:(this.iterations*this.innerIterations) if this.console_output == CONSOLE_OUTPUT_VERBOSE && i > 1
if this.console_output == CONSOLE_OUTPUT_VERBOSE && i > 1 printConcentrationsConsole(this)
printConcentrationsConsole(this)
end
if this.csv_output == CSV_OUTPUT_VERBOSE || this.csv_output == CSV_OUTPUT_XTREME
printConcentrationsCSV(this, filename)
end
# FTCS!(this.grid, this.bc, this.timestep)
end end
elseif this.approach == BTCS_APPROACH
for i = 1:this.iterations
if this.console_output == CONSOLE_OUTPUT_VERBOSE && i > 1
printConcentrationsConsole(this)
end
if this.csv_output == CSV_OUTPUT_VERBOSE || this.csv_output == CSV_OUTPUT_XTREME
printConcentrationsCSV(this, filename)
end
BTCS!(this.grid, this.bc, this.timestep, this.solver) if this.csv_output == CSV_OUTPUT_VERBOSE || this.csv_output == CSV_OUTPUT_XTREME
printConcentrationsCSV(this, filename)
end end
elseif this.approach == CRANK_NICOLSON_APPROACH
beta::T = 0.5
# TODO this implementation is very inefficient! BTCS!(this.grid, this.bc, this.timestep, this.solver)
# a separate implementation that sets up a specific tridiagonal matrix
# for Crank-Nicolson would be better
for i = 1:(this.iterations*this.innerIterations)
if this.console_output == CONSOLE_OUTPUT_VERBOSE && i > 1
printConcentrationsConsole(this)
end
if this.csv_output == CSV_OUTPUT_VERBOSE || this.csv_output == CSV_OUTPUT_XTREME
printConcentrationsCSV(this, filename)
end
concentrations::AbstractMatrix{T} = getGridConcentrations(this.grid)
# FTCS!(this.grid, this.bc, this.timestep)
concentrationsFTCS::AbstractMatrix{T} = getGridConcentrations(this.grid)
setGridConcentrations!(this.grid, concentrations)
BTCS!(this.grid, this.bc, this.timestep, THOMAS_ALGORITHM_SOLVER)
concentrationsResult::AbstractMatrix{T} =
beta * concentrationsFTCS + (1 - beta) * getGridConcentrations(this.grid)
setGridConcentrations!(this.grid, concentrationsResult)
end
end end
if this.console_output != CONSOLE_OUTPUT_OFF if this.console_output != CONSOLE_OUTPUT_OFF
...@@ -374,8 +256,7 @@ function runSimulation!(this::Simulation{T}) where {T} ...@@ -374,8 +256,7 @@ function runSimulation!(this::Simulation{T}) where {T}
end end
if this.time_measure != TIME_MEASURE_OFF if this.time_measure != TIME_MEASURE_OFF
approachString::String = this.approach_names[this.approach]
dimString::String = @sprintf "%dD" getGridDim(this.grid) dimString::String = @sprintf "%dD" getGridDim(this.grid)
println(stdout, approachString, dimString, ":: run() finished in ", Int(milliseconds), "ms") println(stdout, dimString, ":: run() finished in ", Int(milliseconds), "ms")
end end
end end
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment