Skip to content

Commit

Permalink
Merge pull request #4 from BattMoTeam/p2d-perf-moyner
Browse files Browse the repository at this point in the history 
Improve performance of P2D models
  • Loading branch information
@moyner
moyner authored Jun 27, 2023
2 parents 3b3ecaa + 530457d commit 776cd0f
Show file tree
Hide file tree
Showing 6 changed files with 82 additions and 67 deletions.
4 changes: 3 additions & 1 deletion Project.toml
View file
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "BattMo"
uuid = "6f0c0536-3c2c-4762-a987-c605a8a6f898"
authors = ["Olav Møyner <olav.moyner@gmail.com>"]
version = "0.1.1"
version = "0.1.2"

[deps]
JLD2 = "033835bb-8acc-5ee8-8aae-3f567f8a3819"
Expand All @@ -10,6 +10,7 @@ Jutul = "2b460a1a-8a2b-45b2-b125-b5c536396eb9"
MAT = "23992714-dd62-5051-b70f-ba57cb901cac"
Polynomials = "f27b6e38-b328-58d1-80ce-0feddd5e7a45"
PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a"
StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
Tullio = "bc48ee85-29a4-5162-ae0b-a64e1601d4bc"

Expand All @@ -20,5 +21,6 @@ Jutul = "0.2.10"
MAT = "0.10.4"
Polynomials = "3.2.11"
PrecompileTools = "1.1.1"
StaticArrays = "1.5.17"
Tullio = "0.3.5"
julia = "1.8"
2 changes: 2 additions & 0 deletions src/BattMo.jl
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,7 @@
module BattMo
using PrecompileTools
using StaticArrays

import JSON
import Jutul:
number_of_cells, number_of_faces,
Expand Down
89 changes: 48 additions & 41 deletions src/models/activematerial.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -5,8 +5,8 @@ const ActiveMaterialParameters = JutulStorage

abstract type SolidDiffusionDiscretization end

struct P2Ddiscretization <: SolidDiffusionDiscretization
data::Dict{Symbol, Any}
struct P2Ddiscretization{T} <: SolidDiffusionDiscretization
data::T
# At the moment the following keys are included :
# N::Integer # Discretization size for solid diffusion
# rp::Real # Particle radius
Expand All @@ -18,8 +18,8 @@ end

struct NoParticleDiffusion <: SolidDiffusionDiscretization end

struct ActiveMaterial{D} <: ElectroChemicalComponent where {D<:SolidDiffusionDiscretization}
params::ActiveMaterialParameters
struct ActiveMaterial{D, T} <: ElectroChemicalComponent where {D<:SolidDiffusionDiscretization, T<:ActiveMaterialParameters}
params::T
# At the moment the following keys are include
# - ocp_func::F where {F <: Function}
# - n_charge_carriers
Expand All @@ -30,6 +30,9 @@ struct ActiveMaterial{D} <: ElectroChemicalComponent where {D<:SolidDiffusionDis
discretization::D
end

const ActiveMaterialP2D{D, T} = ActiveMaterial{D, T}
const ActiveMaterialNoParticleDiffusion{T} = ActiveMaterial{NoParticleDiffusion, T}

struct Ocp <: ScalarVariable end
struct DiffusionCoef <: ScalarVariable end
struct ReactionRateConst <: ScalarVariable end
Expand All @@ -46,16 +49,18 @@ Jutul.local_discretization(::SolidDiffusionBc, i) = nothing
const ActiveMaterialModel = SimulationModel{O, S} where {O<:JutulDomain, S<:ActiveMaterial}

## Create ActiveMaterial with full p2d solid diffusion
function ActiveMaterial{P2Ddiscretization}(params::ActiveMaterialParameters, rp, N, D)
function ActiveMaterialP2D(params::ActiveMaterialParameters, rp, N, D)
data = setupSolidDiffusionDiscretization(rp, N, D)
discretization = P2Ddiscretization(data)
return ActiveMaterial{P2Ddiscretization}(params, discretization)
params = Jutul.convert_to_immutable_storage(params)
return ActiveMaterial{typeof(discretization), typeof(params)}(params, discretization)
end

## Create ActiveMaterial with no solid diffusion
function ActiveMaterial{NoParticleDiffusion}(params::ActiveMaterialParameters)
function ActiveMaterialNoParticleDiffusion(params::ActiveMaterialParameters)
discretization = NoParticleDiffusion()
return ActiveMaterial{NoParticleDiffusion}(params, discretization)
params = Jutul.convert_to_immutable_storage(params)
return ActiveMaterial{NoParticleDiffusion, typeof(params)}(params, discretization)
end


Expand All @@ -68,19 +73,19 @@ function Base.getindex(disc::P2Ddiscretization, key::Symbol)
return disc.data[key]
end

function discretisation_type(system::ActiveMaterial{P2Ddiscretization})
function discretisation_type(system::ActiveMaterialP2D)
return :P2Ddiscretization
end

function discretisation_type(system::ActiveMaterial{NoParticleDiffusion})
function discretisation_type(system::ActiveMaterialNoParticleDiffusion)
return :NoParticleDiffusion
end

function discretisation_type(model::ActiveMaterialModel)
discretisation_type(model.system)
end

function solid_diffusion_discretization_number(system::ActiveMaterial{P2Ddiscretization})
function solid_diffusion_discretization_number(system::ActiveMaterialP2D)
return system.discretization[:N]
end

Expand Down Expand Up @@ -108,7 +113,7 @@ function setupSolidDiffusionDiscretization(rp, N, D)
for i = 1 : N + 1
hT[i] = (4*pi*rf[i]^2/(dr/2))
end

div = Vector{Tuple{Int64, Int64, Float64}}(undef, 2*(N - 1))

k = 1
Expand All @@ -118,25 +123,27 @@ function setupSolidDiffusionDiscretization(rp, N, D)
div[k] = (j + 1, j, -1)
k += 1
end

hT = SVector(hT...)
div = SVector(div...)
vols = SVector(vols...)

data = Dict(:N => N ,
:D => D ,
:rp => rp ,
:hT => hT ,
:vols => vols,
:div => div ,
)

return data

# Convert to concrete type
return NamedTuple(pairs(data))
end

#####################################################
## We setup the case with full P2d discretization
####################################################

function select_primary_variables!(S,
system::ActiveMaterial{P2Ddiscretization},
system::ActiveMaterialP2D,
model::SimulationModel
)
S[:Phi] = Phi()
Expand All @@ -156,7 +163,7 @@ function degrees_of_freedom_per_entity(model::ActiveMaterialModel,
end

function select_parameters!(S,
system::ActiveMaterial{P2Ddiscretization},
system::ActiveMaterialP2D,
model::SimulationModel)

S[:Temperature] = Temperature()
Expand All @@ -166,7 +173,7 @@ function select_parameters!(S,
end

function select_secondary_variables!(S,
system::ActiveMaterial{P2Ddiscretization},
system::ActiveMaterialP2D,
model::SimulationModel
)
S[:Charge] = Charge()
Expand All @@ -177,7 +184,7 @@ function select_secondary_variables!(S,
end

function select_equations!(eqs,
system::ActiveMaterial{P2Ddiscretization},
system::ActiveMaterialP2D,
model::SimulationModel
)

Expand All @@ -197,7 +204,7 @@ function Jutul.number_of_equations_per_entity(model::ActiveMaterialModel, ::Soli
end

function select_minimum_output_variables!(out ,
system::ActiveMaterial{P2Ddiscretization},
system::ActiveMaterialP2D,
model::SimulationModel
)
push!(out, :Charge)
Expand All @@ -209,10 +216,10 @@ end
@jutul_secondary(
function update_vocp!(Ocp,
tv::Ocp,
model:: SimulationModel{<:Any, ActiveMaterial{P2Ddiscretization}, <:Any, <:Any},
model:: SimulationModel{<:Any, ActiveMaterialP2D{D, T}, <:Any, <:Any},
Cs,
ix
)
) where {D, T}
ocp_func = model.system.params[:ocp_func]
cmax = model.system.params[:maximum_concentration]
refT = 298.15
Expand All @@ -226,10 +233,10 @@ end
@jutul_secondary(
function update_reaction_rate!(ReactionRateConst ,
tv::ReactionRateConst ,
model::SimulationModel{<:Any, ActiveMaterial{P2Ddiscretization}, <:Any, <:Any},
model::SimulationModel{<:Any, ActiveMaterialP2D{D, T}, <:Any, <:Any},
Cs ,
ix
)
) where {D, T}
rate_func = model.system.params[:reaction_rate_constant_func]
refT = 298.15
for cell in ix
Expand All @@ -241,9 +248,9 @@ end
@jutul_secondary(
function update_solid_diffusion_flux!(SolidDiffFlux,
tv::SolidDiffFlux,
model::SimulationModel{<:Any, ActiveMaterial{P2Ddiscretization}, <:Any, <:Any},
model::SimulationModel{<:Any, ActiveMaterialP2D{D, T}, <:Any, <:Any},
Cp,
ix)
ix) where {D, T}
s = model.system
for cell in ix
@inbounds @views update_solid_flux!(SolidDiffFlux[:, cell], Cp[:, cell], s)
Expand All @@ -252,7 +259,7 @@ end
)


function update_solid_flux!(flux, Cp, system::ActiveMaterial{P2Ddiscretization})
function update_solid_flux!(flux, Cp, system::ActiveMaterialP2D)
# compute lithium flux in particle, using harmonic averaging. At the moment D has a constant value within particle
# but this is going to change.

Expand Down Expand Up @@ -286,15 +293,15 @@ function Jutul.update_equation_in_entity!(eq_buf ,
ldisc = Nothing)

disc = model.system.discretization
N = disc[:N]
N = length(eq_buf)
hT = disc[:hT]
vols = disc[:vols]
div = disc[:div]
D = disc[:D]

Cp = state.Cp[:, self_cell]
Cp0 = state0.Cp[:, self_cell]
flux = state.SolidDiffFlux[:, self_cell]
Cp = @views state.Cp[:, self_cell]
Cp0 = @views state0.Cp[:, self_cell]
flux = @views state.SolidDiffFlux[:, self_cell]
Cs = state.Cs[self_cell]

for i = 1 : N
Expand Down Expand Up @@ -339,7 +346,7 @@ end


function select_primary_variables!(S,
system::ActiveMaterial{NoParticleDiffusion},
system::ActiveMaterialNoParticleDiffusion,
model::SimulationModel
)
S[:Phi] = Phi()
Expand All @@ -348,7 +355,7 @@ function select_primary_variables!(S,
end

function select_secondary_variables!(S,
system::ActiveMaterial{NoParticleDiffusion},
system::ActiveMaterialNoParticleDiffusion,
model::SimulationModel
)

Expand All @@ -360,7 +367,7 @@ function select_secondary_variables!(S,
end

function select_parameters!(S,
system::ActiveMaterial{NoParticleDiffusion},
system::ActiveMaterialNoParticleDiffusion,
model::SimulationModel)

S[:Temperature] = Temperature()
Expand All @@ -370,7 +377,7 @@ function select_parameters!(S,
end

function select_equations!(eqs,
system::ActiveMaterial{NoParticleDiffusion},
system::ActiveMaterialNoParticleDiffusion,
model::SimulationModel
)

Expand All @@ -382,7 +389,7 @@ end


function select_minimum_output_variables!(out ,
system::ActiveMaterial{NoParticleDiffusion},
system::ActiveMaterialNoParticleDiffusion,
model::SimulationModel
)
push!(out, :Charge)
Expand All @@ -397,10 +404,10 @@ end
@jutul_secondary(
function update_vocp!(Ocp ,
tv::Ocp ,
model::M,
model::SimulationModel{<:Any, ActiveMaterialNoParticleDiffusion{T}, <:Any, <:Any},
C ,
ix
) where {M <: SimulationModel{<:Any, ActiveMaterial{NoParticleDiffusion}, <:Any, <:Any}}
) where T

ocp_func = model.system.params[:ocp_func]
cmax = model.system.params[:maximum_concentration]
Expand All @@ -417,10 +424,10 @@ end
@jutul_secondary(
function update_reaction_rate!(ReactionRateConst ,
tv::ReactionRateConst,
model::SimulationModel{<:Any, ActiveMaterial{NoParticleDiffusion}, <:Any, <:Any},
model::SimulationModel{<:Any, ActiveMaterialNoParticleDiffusion{T}, <:Any, <:Any},
C ,
ix
)
) where T
rate_func = model.system.params[:reaction_rate_constant_func]
refT = 298.15
for i in ix
Expand Down
6 changes: 2 additions & 4 deletions src/models/battery_cross_terms.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -87,13 +87,13 @@ function reaction_rate(phi_a ,
n = activematerial.params[:n_charge_carriers]
cmax = activematerial.params[:maximum_concentration]
vsa = activematerial.params[:volumetric_surface_area]

eta = (phi_a - phi_e - ocp)
th = 1e-3*cmax
j0 = R0*regularized_sqrt(c_e*(cmax - c_a)*c_a, th)*n*FARADAY_CONST
R = vsa*butler_volmer_equation(j0, 0.5, n, eta, T)

return R./(n*FARADAY_CONST)
return R/(n*FARADAY_CONST)

end

Expand Down Expand Up @@ -159,9 +159,7 @@ function Jutul.update_cross_term_in_entity!(out ,
@assert cs == :Charge
v = 1.0*vols*R*n*FARADAY_CONST
end

out[] = -v

end


Expand Down
10 changes: 4 additions & 6 deletions src/models/elyte.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -154,7 +154,7 @@ end
t = transference(sys)
F = FARADAY_CONST
for i in ix
@inbounds chemCoef[i] = 1/F*(1 - t)*Conductivity[i]*2*DmuDc[i]
@inbounds chemCoef[i] = 1.0/F*(1.0 - t)*Conductivity[i]*2.0*DmuDc[i]
end
end

Expand All @@ -172,10 +172,10 @@ function Jutul.face_flux!(::T, c, other, face, face_sign, eq::ConservationLaw{:C
end


function Jutul.face_flux!(::T, c, other, face, face_sign, eq::ConservationLaw{:Mass, <:Any}, state, model::ElectrolyteModel, dt, flow_disc) where T
function Jutul.face_flux!(q::T, c, other, face, face_sign, eq::ConservationLaw{:Mass, <:Any}, state, model::ElectrolyteModel, dt, flow_disc) where T

t = transference(model.system)
z = 1
z = 1.0
F = FARADAY_CONST

@inbounds trans = state.ECTransmissibilities[face]
Expand All @@ -187,7 +187,5 @@ function Jutul.face_flux!(::T, c, other, face, face_sign, eq::ConservationLaw{:M
j = j - jchem*(1.0)

massFlux = diffFlux + t/(z*F)*j

return T(massFlux)

return setindex(q, massFlux, 1)::T
end
Loading

2 comments on commit 776cd0f

@moyner
Copy link
Member Author

@moyner moyner commented on 776cd0f Jun 27, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@JuliaRegistrator register

@JuliaRegistrator
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Registration pull request created: JuliaRegistries/General/86343

After the above pull request is merged, it is recommended that a tag is created on this repository for the registered package version.

This will be done automatically if the Julia TagBot GitHub Action is installed, or can be done manually through the github interface, or via:

git tag -a v0.1.2 -m "<description of version>" 776cd0f83aebf85eff80f1e2c786af7d1a9b7e5a
git push origin v0.1.2

Please sign in to comment.