@@ -20,7 +20,7 @@ fn = string(dirname(pathof(BattMo)), "/../test/data/jsonfiles/", name, ".json")
2020inputparams = readBattMoJsonInputFile (fn)
2121
2222
23- simple = true
23+ simple = false
2424if (! simple)
2525
2626 facx = 1
@@ -30,29 +30,28 @@ if(!simple)
3030
3131 fn = string (dirname (pathof (BattMo)), " /../test/data/jsonfiles/3d_demo_geometry.json" )
3232 inputparams_geometry = readBattMoJsonInputFile (fn)
33- inputparams_geometry. dict[" Geometry" ][" Nh" ] *= facy
34- inputparams_geometry. dict[" Geometry" ][" Nw" ] *= facx
35- inputparams_geometry. dict[" Separator" ][" N" ] *= facz
36- inputparams_geometry. dict[" PositiveElectrode" ][" Coating" ][" N" ] *= facz
37- inputparams_geometry. dict[" NegativeElectrode" ][" Coating" ][" N" ] *= facz
33+ inputparams_geometry. dict[" Geometry" ][" Nh" ] *= facy
34+ inputparams_geometry. dict[" Geometry" ][" Nw" ] *= facx
35+ inputparams_geometry. dict[" Separator" ][" N" ] *= facz
36+ inputparams_geometry. dict[" PositiveElectrode" ][" Coating" ][" N" ] *= facz
37+ inputparams_geometry. dict[" NegativeElectrode" ][" Coating" ][" N" ] *= facz
3838 inputparams_geometry. dict[" PositiveElectrode" ][" CurrentCollector" ][" tab" ][" Nh" ] *= facy
3939 inputparams_geometry. dict[" PositiveElectrode" ][" CurrentCollector" ][" tab" ][" Nw" ] *= facx
4040 inputparams_geometry. dict[" NegativeElectrode" ][" CurrentCollector" ][" tab" ][" Nh" ] *= facy
4141 inputparams_geometry. dict[" NegativeElectrode" ][" CurrentCollector" ][" tab" ][" Nw" ] *= facx
42- inputparams_geometry. dict[" NegativeElectrode" ][" CurrentCollector" ][" N" ] *= facz
43- inputparams_geometry. dict[" PositiveElectrode" ][" CurrentCollector" ][" N" ] *= facz
42+ inputparams_geometry. dict[" NegativeElectrode" ][" CurrentCollector" ][" N" ] *= facz
43+ inputparams_geometry. dict[" PositiveElectrode" ][" CurrentCollector" ][" N" ] *= facz
44+
4445else
4546 # fn = string(dirname(pathof(BattMo)), "/../test/data/jsonfiles/3d_demo_geometry.json")
4647 fn = string (dirname (pathof (BattMo)), " /../test/data/jsonfiles/1D_geometry.json" )
4748 inputparams_geometry_org = readBattMoJsonInputFile (fn)
49+ inputparams_geometry = deepcopy (inputparams_geometry_org)
50+ inputparams_geometry. dict[" include_current_collectors" ] = false
4851end
4952
50- inputparams_geometry = deepcopy (inputparams_geometry_org)
51-
52- inputparams_geometry. dict[" include_current_collectors" ] = false
5353inputparams = mergeInputParams (deepcopy (inputparams_geometry), inputparams)
5454
55-
5655# ###########################
5756# setup and run simulation #
5857# ###########################
@@ -82,18 +81,28 @@ atol = 1e-5*fac # seems important
8281max_it = 100
8382verbose = 10
8483
85- varpreconds = Vector {BattMo.VariablePrecond} ()
84+ # We combine two preconditioners. One working on a subset of variables and equations (we call it block-preconditioner)
85+ # and the other for the full system
8686
87+ # We first setup the block preconditioners. They are given as a list and applied separatly. Preferably, they
88+ # should be orthogonal
89+ varpreconds = Vector {BattMo.VariablePrecond} ()
8790push! (varpreconds, BattMo. VariablePrecond (Jutul. AMGPreconditioner (:ruge_stuben ),:Phi ,:charge_conservation , nothing ))
8891# push!(varpreconds,BattMo.VariablePrecond(Jutul.ILUZeroPreconditioner(),:Cp,:mass_conservation, [:PeAm,:NeAm]))
8992# push!(varpreconds,BattMo.VariablePrecond(Jutul.AMGPreconditioner(:ruge_stuben),:C,:mass_conservation, [:Elyte]))
90- g_varprecond = BattMo. VariablePrecond (Jutul. ILUZeroPreconditioner (),:Global ,:Global ,nothing )
93+
94+ # We setup the global preconditioner
95+ g_varprecond = BattMo. VariablePrecond (Jutul. ILUZeroPreconditioner (), :Global , :Global ,nothing )
9196
9297params = Dict ()
98+ # Type of method used for the block preconditioners. Here "block" means separatly (other options can be found
99+ # BatteryGeneralPreconditione)
93100params[" method" ] = " block"
101+ # Option for post- and pre-solve of the control system.
94102params[" post_solve_control" ] = true
95103params[" pre_solve_control" ] = true
96104
105+ # We setup the preconditioner, which combines both the block and global preconditioners
97106prec = BattMo. BatteryGeneralPreconditioner (varpreconds, g_varprecond, params)
98107# prec = Jutul.ILUZeroPreconditioner()
99108
@@ -102,7 +111,7 @@ cfg[:linear_solver] = GenericKrylov(solver, verbose = verbose,
102111 relative_tolerance = rtol,
103112 absolute_tolerance = atol* 1e-20 ,# # may skip linear iterations all to getter.
104113 max_iterations = max_it)
105- cfg[:extra_timing ] = true
114+ cfg[:extra_timing ] = true
106115
107116# Perform simulation
108117states, reports = simulate (state0, simulator, timesteps; forces = forces, config = cfg)
@@ -162,101 +171,105 @@ display(f)
162171# ###########################################
163172# plot potential on grid at last time step #
164173# ###########################################
165- do_plot = false
174+
175+ do_plot = true
176+
166177if (do_plot)
167- state = states[10 ]
168178
169- setups = ((:PeCc , :PeAm , " positive" ),
170- (:NeCc , :NeAm , " negative" ))
179+ state = states[10 ]
171180
181+ setups = ((:PeCc , :PeAm , " positive" ),
182+ (:NeCc , :NeAm , " negative" ))
172183
173- for setup in setups
174184
175- f3D = Figure (size = (600 , 650 ))
176- ax3d = Axis3 (f3D[1 , 1 ];
177- title = " Potential in $(setup[3 ]) electrode (coating and active material)" )
185+ for setup in setups
178186
179- am = setup[1 ]
180- cc = setup[2 ]
181-
182- maxPhi = maximum ([maximum (state[cc][:Phi ]), maximum (state[am][:Phi ])])
183- minPhi = minimum ([minimum (state[cc][:Phi ]), minimum (state[am][:Phi ])])
187+ f3D = Figure (size = (600 , 650 ))
188+ ax3d = Axis3 (f3D[1 , 1 ];
189+ title = " Potential in $(setup[3 ]) electrode (coating and active material)" )
184190
185- colorrange = [0 , maxPhi - minPhi]
191+ am = setup[1 ]
192+ cc = setup[2 ]
193+
194+ maxPhi = maximum ([maximum (state[cc][:Phi ]), maximum (state[am][:Phi ])])
195+ minPhi = minimum ([minimum (state[cc][:Phi ]), minimum (state[am][:Phi ])])
186196
187- components = [am, cc]
188- for component in components
189- g = model[component]. domain. representation
190- phi = state[component][:Phi ]
191- Jutul. plot_cell_data! (ax3d, g, phi .- minPhi; colormap = :viridis , colorrange = colorrange)
192- end
197+ colorrange = [0 , maxPhi - minPhi]
193198
194- cbar = GLMakie. Colorbar (f3D[1 , 2 ];
195- colormap = :viridis ,
196- colorrange = colorrange .+ minPhi,
197- label = " potential" )
198- display (GLMakie. Screen (), f3D)
199+ components = [am, cc]
200+ for component in components
201+ g = model[component]. domain. representation
202+ phi = state[component][:Phi ]
203+ Jutul. plot_cell_data! (ax3d, g, phi .- minPhi; colormap = :viridis , colorrange = colorrange)
204+ end
199205
200- end
206+ cbar = GLMakie. Colorbar (f3D[1 , 2 ];
207+ colormap = :viridis ,
208+ colorrange = colorrange .+ minPhi,
209+ label = " potential" )
210+ display (GLMakie. Screen (), f3D)
201211
202- setups = ((:PeAm , " positive" ),
203- (:NeAm , " negative" ))
212+ end
204213
205- for setup in setups
214+ setups = ((:PeAm , " positive" ),
215+ (:NeAm , " negative" ))
206216
207- f3D = Figure (size = (600 , 650 ))
208- ax3d = Axis3 (f3D[1 , 1 ];
209- title = " Surface concentration in $(setup[2 ]) electrode" )
217+ for setup in setups
210218
211- component = setup[1 ]
212-
213- cs = state[component][:Cs ]
214- maxcs = maximum (cs)
215- mincs = minimum (cs)
219+ f3D = Figure (size = (600 , 650 ))
220+ ax3d = Axis3 (f3D[1 , 1 ];
221+ title = " Surface concentration in $(setup[2 ]) electrode" )
216222
217- colorrange = [0 , maxcs - mincs]
223+ component = setup[1 ]
224+
225+ cs = state[component][:Cs ]
226+ maxcs = maximum (cs)
227+ mincs = minimum (cs)
218228
219- g = model[component]. domain. representation
220- Jutul. plot_cell_data! (ax3d, g, cs .- mincs;
221- colormap = :viridis ,
222- colorrange = colorrange)
229+ colorrange = [0 , maxcs - mincs]
223230
224- cbar = GLMakie. Colorbar (f3D[1 , 2 ];
225- colormap = :viridis ,
226- colorrange = colorrange .+ mincs,
227- label = " concentration" )
228- display (GLMakie. Screen (), f3D)
231+ g = model[component]. domain. representation
232+ Jutul. plot_cell_data! (ax3d, g, cs .- mincs;
233+ colormap = :viridis ,
234+ colorrange = colorrange)
229235
230- end
236+ cbar = GLMakie. Colorbar (f3D[1 , 2 ];
237+ colormap = :viridis ,
238+ colorrange = colorrange .+ mincs,
239+ label = " concentration" )
240+ display (GLMakie. Screen (), f3D)
231241
242+ end
232243
233- setups = ((:C , " concentration" ),
234- (:Phi , " potential" ))
235244
236- for setup in setups
245+ setups = ((:C , " concentration" ),
246+ (:Phi , " potential" ))
237247
238- f3D = Figure (size = (600 , 650 ))
239- ax3d = Axis3 (f3D[1 , 1 ];
240- title = " $(setup[2 ]) in electrolyte" )
248+ for setup in setups
241249
242- var = setup[1 ]
243-
244- val = state[:Elyte ][var]
245- maxval = maximum (val)
246- minval = minimum (val)
250+ f3D = Figure (size = (600 , 650 ))
251+ ax3d = Axis3 (f3D[1 , 1 ];
252+ title = " $(setup[2 ]) in electrolyte" )
247253
248- colorrange = [0 , maxval - minval]
254+ var = setup[1 ]
255+
256+ val = state[:Elyte ][var]
257+ maxval = maximum (val)
258+ minval = minimum (val)
249259
250- g = model[:Elyte ]. domain. representation
251- Jutul. plot_cell_data! (ax3d, g, val .- minval;
252- colormap = :viridis ,
253- colorrange = colorrange)
260+ colorrange = [0 , maxval - minval]
254261
255- cbar = GLMakie. Colorbar (f3D[1 , 2 ];
256- colormap = :viridis ,
257- colorrange = colorrange .+ minval,
258- label = " $(setup[2 ]) " )
259- display (GLMakie. Screen (), f3D)
262+ g = model[:Elyte ]. domain. representation
263+ Jutul. plot_cell_data! (ax3d, g, val .- minval;
264+ colormap = :viridis ,
265+ colorrange = colorrange)
260266
261- end
267+ cbar = GLMakie. Colorbar (f3D[1 , 2 ];
268+ colormap = :viridis ,
269+ colorrange = colorrange .+ minval,
270+ label = " $(setup[2 ]) " )
271+ display (GLMakie. Screen (), f3D)
272+
273+ end
274+
262275end
0 commit comments