1 function model = advection_fv_output_model_hp(params)
2 %
function model = advection_fv_output_model(params);
4 % model
for the
new advection model with functional output
5 % two overlapping velocity fields and boundary value modification
6 % are the parameters. Output functional is average concentration
7 % in subdomain. Pure
explicit discretization of convection with ldg
10 % example of a model not as
class but as structure plus function pointers
11 % goal: library works with structure- or
class-models
13 % possible fields of params:
14 % coarse_factor: coarsening factor between 1 and 8
15 % time-steps, gridsize are scaled down with this.
17 % (To be used as base model
for initial results of RBEvolOpt...)
19 % B. Haasdonk 26.8.2009
22 if (nargin ==1) & (isfield(params,
'coarse_factor'))
23 coarse_factor = params.coarse_factor;
28 % specification of the time information
30 model.nt = 2048/coarse_factor;
31 %model.save_time_indices = 0:16/coarse_factor:model.nt;
32 model.save_time_indices = 0:8/coarse_factor:model.nt;
33 disp(
'nt to be adjusted later!');
34 model.dt = model.T/model.nt;
39 % disp(
'model with debugging turned on.');
42 model.error_estimation = 1;
45 %model.grid_initfile =
'rectangle_triagrid.mat';
46 model.gridtype =
'triagrid';
47 model.xnumintervals = 256/coarse_factor;
48 model.ynumintervals = 128/coarse_factor;
51 % set completely dirichlet boundary
52 model.bnd_rect_corner1 = [0,0]-eps;
53 model.bnd_rect_corner2 = [2,1]+eps;
54 %model.opts.bnd_rect_index = [-1];
55 model.element_quadrature = @triaquadrature;
57 model.dim_U = model.xnumintervals * model.ynumintervals * 2;
59 % Main global
function pointers expected in every model
60 model.gen_model_data = @lin_evol_gen_model_data;
61 model.gen_detailed_data = @lin_evol_gen_detailed_data;
62 model.gen_reduced_data = @lin_evol_gen_reduced_data;
63 model.detailed_simulation = @lin_evol_detailed_simulation;
64 model.rb_simulation = @lin_evol_rb_simulation;
66 model.plot_sim_data = @lin_evol_plot_sim_data;
69 model.distance_function = @euclidean_distance;
70 model.RB_extension_algorithm_hp = @RB_extension_PCA_fixspace_flexible_hp;
71 model.PCA_fixspace_evalues = @PCA_fixspace_evalues;
74 % Dirichlet and Initial data
75 model.dirichlet_values_ptr = @dirichlet_values_affine_decomposed;
76 model.dirichlet_values_coefficients_ptr = @my_dirichlet_values_coefficients;
77 model.dirichlet_values_components_ptr = @my_dirichlet_values_components;
79 %model.dirichlet_values_coefficients_ptr = @(params) 1;
80 %model.dirichlet_values_components_ptr = @(glob,params) {ones(1, ...
82 model.rb_operators = @lin_evol_rb_operators;
83 model.error_norm =
'l2';
84 model.L_I_inv_norm_bound = 1;
85 model.L_E_norm_bound = 1;
86 model.rb_init_data_basis = @RB_init_data_basis;
88 model.init_values_ptr = @init_values_affine_decomposed;
89 model.init_values_coefficients_ptr = @my_dirichlet_values_coefficients_t0;
90 model.init_values_components_ptr = @my_dirichlet_values_components;
91 model.init_values_algorithm = @fv_init_values;
92 model.cone_number = 3; % number of components = cones
93 model.cone_range = [0,1]; % x-range of cone-support
94 model.cone_weight = 1; % leftmost cone with full weight
95 model.velocity_ptr = @velocity_affine_decomposed;
96 model.velocity_coefficients_ptr = @my_velocity_coefficients;
97 model.velocity_components_ptr = @my_velocity_components;
98 %model.vx_weight = 1.0;
99 model.vx_weight = 0.75;
100 %model.vy_weight = 1.0;
101 model.vy_weight = 0.75;
103 model.set_rb_in_detailed_data = @(detailed_data,RB) ...
104 setfield(detailed_data,
'RB',RB);
105 model.get_rb_size = @(model,detailed_data)size(detailed_data.RB,2);
106 % set parameter setting
function
107 model.set_time = @set_time_default; % not set_time!!!
108 model.set_mu = @set_mu_default;
109 model.get_mu = @get_mu_default;
111 % perhaps these are redundant later...
112 model.divclean_mode = 0;
113 model.flux_quad_degree = 1;
114 model.flux_linear = 1;
116 %model.init_values_algorithm = @disc_init_values;
117 model.init_values_qdeg = 0;
118 model.pdeg = 0; % FV schemes with piecewise constant ansatz functions
119 model.evaluate_basis = @fv_evaluate_basis;
120 model.l2project = @l2project;
121 model.local_mass_matrix = @fv_local_mass_matrix_tria; % triangular grid
122 model.mass_matrix = @fv_mass_matrix;
123 model.ndofs_per_element = 1;
126 %
for use of old datafunctions, use the following:
127 %model.operators_ptr = @fv_operators_implicit_explicit;
128 %model.operators_conv_explicit = @fv_operators_conv_explicit;
129 %model.operators_conv_implicit = @fv_operators_conv_implicit;
130 %model.operators_diff_explicit = @fv_operators_diff_explicit;
131 %model.operators_diff_implicit = @fv_operators_diff_implicit;
132 %model.operators_neumann_explicit = @fv_operators_neumann_explicit_old;
133 %model.operators_neumann_implicit = @fv_operators_neumann_implicit;
135 model.operators_ptr = @fv_operators_implicit_explicit;
136 %
for use of
new data
function access, use now
138 model.operators_conv_implicit = @fv_operators_zero;
139 model.operators_diff_explicit = @fv_operators_zero;
140 model.operators_diff_implicit = @fv_operators_zero;
141 model.operators_neumann_implicit = @fv_operators_zero;
142 model.operators_neumann_explicit = @fv_operators_zero;
144 model.orthonormalize= @model_orthonormalize_gram_schmidt; %
new
148 %model.flux_linear = 1;
149 model.data_const_in_time = 0; % time varying data...
150 %model.diffusivity_ptr = ...;
151 %model.neumann_value_ptr = @...
152 model.affinely_decomposed = 1; % data functions allow affine
153 % parameter decomposition
155 model.compute_output_functional = 1; % turn on computation of output
156 model.output_function_ptr = @output_function_box_mean;
157 %model.sbox_xmin = 1;
161 %model.sbox_ymax = 0.5;
162 model.sbox_ymax = 0.5;
167 %model.mu_names = {
'cone_weight',
'vx_weight',
'vy_weight'};
168 model.mu_names = {
'vx_weight',
'vy_weight'};
169 %model.mu_ranges = {[0 1],[0 1],[0,1]};
170 model.mu_ranges = {[0 1],[0 1]};%,[0 1]};
173 %model.mu_names = {
'cone_weight',
'vx_weight',
'vy_weight'};
174 %model.mu_ranges = {[0 1],[0 1],[0,1]};
176 model.random_seed=1234;
177 %model.RB_numintervals = [2,2,2];
178 %model.RB_numintervals = [2];
179 model.RB_generation_mode =
'greedy_refined';%
'greedy_uniform_fixed'; %
180 model.RB_refinement_mode =
'adaptive';
181 model.RB_error_indicator =
'estimator';%
'error';%
'estimator';
182 model.RB_detailed_train_savepath = fullfile(rbmatlabtemp,
'dump');
184 model.RB_stop_epsilon = 1e-2;%1e-1;
185 model.RB_stop_timeout = 24*60*60; % 1 minute
186 model.RB_stop_Nmax = 20;
187 model.nr_extension_modes = 1;
188 %fields
for adaptive training set extension:
189 %model.RB_stop_max_val_train_ratio = 1;
190 %model.RB_refinement_theta = 0.2;
191 %model.RB_val_rand_seed = 543;
192 %model.RB_M_val_size = 10;
193 %model.RB_max_refinement_level = 5;
194 %model.RB_element_indicator_mode=
'nodes_cogs_skippedrefs';
195 %model.RB_element_indicator_s_max = 1;
196 %model.RB_mixed_sensitivity_basis_generation =0;
198 model.get_inner_product_matrix = @(model_data)model_data.W;
199 model.inner_product_matrix_algorithm = @(model,model_data)model_data.W;
200 model.get_estimator_from_sim_data = @(sim_data)sim_data.Delta(end);
201 model.get_dofs_from_sim_data = @(sim_data) sim_data.U;
202 model.filecache_ignore_fields_in_model = {
'N',
'Nmax',
'mu_ranges'};
203 model.filecache_ignore_fields_in_detailed_data = {
'RB_info'};
204 model.get_rb_from_detailed_data = @(detailed_data) detailed_data.RB;
205 model.set_rb_in_detailed_data = @(detailed_data,RB) ...
206 setfield(detailed_data,
'RB',RB);
207 model.PCA_fixspace = @PCA_fixspace;
209 %
return; % temporary end
211 %model.L_I_inv_norm_bound = 1; % bounds
for implicit/
explicit operator
212 %model.L_E_norm_bound = 1;
214 %model.neumann_values_ptr = 0;
216 %model.name_init_values =
'decomp_function_ptr';
218 %model.init_values_coefficients_ptr = @my_init_values_coefficients;
219 %model.init_values_components_ptr = @my_init_values_components;
221 %model.name_diffusive_num_flux =
'none';
224 %name_diffusive_num_flux =
'gradient';
225 % precomputed, divcleaned velocity field
226 % name_flux =
'gdl2';
228 % lambda = 2.0e-7; % v = - lambda * grad p
229 % name_convective_num_flux =
'lax-friedrichs';
230 % inner_product_matrix_algorithm = @fv_inner_product_matrix;
233 %model.rb_init_values = @rb_init_values;
235 %% further method pointers, which are specific to model:
236 model.orthonormalize = @model_orthonormalize_gram_schmidt;
237 model.inner_product = @fv_inner_product;
238 %model.PCA = @model_PCA_fixspace;
239 %model.cached_detailed_simulation = @cached_detailed_simulation;
242 %model.use_velocity_matrixfile = 1;
243 %model.divclean_mode =
'none'; % file is already cleaned by optimization
245 % set matrix-file name
for optional generation or reading of file
246 %model.velocity_matrixfile = [
'vel_', model.name_flux,
'_',...
247 % num2str( model.xnumintervals),
'x',...
248 % num2str( model.ynumintervals),...
249 %
'_l',num2str( model.lambda),
'.mat'];%
251 %model.lxf_lambda = 1.0194e+003;
253 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
254 %%%%% auxiliary functions used as pointers above:
255 %%%%% check later,
if they are of general interest to be exported
256 %%%%% as standalone functions
257 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
259 %
new function syntax: global coordinates as rows in glob_coord,
260 % time and parameter variables in params
262 % dirichlet-data: convex combination of equidistant cones, decaying
264 function res = my_dirichlet_values_coefficients_t0(params);
266 res = my_dirichlet_values_coefficients(params);
268 function res = my_dirichlet_values_coefficients(params);
269 % res is a vector of coefficients
270 Q_0 = params.cone_number;
272 max_pos = params.cone_weight * (Q_0-1)+1;
275 res(q) = (1-(max_pos-q))*(1-t) * ((max_pos>=q) && (max_pos < q+1)) ...
276 + (1+(max_pos-q))*(1-t) * ((max_pos>=q-1) && (max_pos < q));
279 function res = my_dirichlet_values_components(glob,params);
280 % res is a cell-array of row-vectors of global evaluations
281 Q_0 = params.cone_number;
283 %delta_cone = 1/(Q_0+1);
284 delta_cone = (params.cone_range(2)-params.cone_range(1))/(Q_0+1);
285 cone_pos_x = delta_cone * (1:Q_0)+ params.cone_range(1);
287 res{q} = 1-min(sqrt((glob(:,1)-cone_pos_x(q)).^2+...
288 (glob(:,2)-1).^2)*(Q_0+1),1);
291 % velocity field: overlap of y- and x parabolic profiles
292 function res = my_velocity_coefficients(params);
293 res = [params.vx_weight, params.vy_weight]*(1-params.t);
295 function res = my_velocity_components(glob,params);
296 resx = [5*(1-glob(:,2).^2),...
297 zeros(size(glob,1),1)];
298 %resy = [zeros(1,size(glob,2));...
299 % -2*((1-glob(1,:).^2 .*(glob(1,:)>=0) & (glob(1,:)<=1)...
301 resy = [zeros(size(glob,1),1),...
302 -1*((4-glob(:,1).^2))];
305 %
function res = my_velocity_coefficients2(params);
308 %
function res = my_velocity_components2(glob,params);
309 %res = {repmat([0;-1],1,size(glob,2))};
function [ sim_data , tictoc ] = load_detailed_simulation(m, savepath, params)
load single trajectory of previously saved results.
function res = intervalquadrature(poldeg, func, varargin)
integration of function func over reference interval == unit interval. by Gaussian quadrature exactly...
function r = verbose(level, message, messageId)
This function displays messages depending on a message-id and/or a level. Aditionally you can set/res...
function [ L_I_diff , bdir_I_diff ] = fv_operators_diff_implicit_gradient(model, model_data, U, NU_ind)
computes diffusion contributions to finite volume time evolution matrices, or their Frechet derivati...
function [ L_E_conv , bdir_E_conv ] = fv_operators_conv_explicit_engquist_osher(model, model_data, U, NU_ind)
computes convection contribution to finite volume time evolution matrices, or their Frechet derivati...
function Udirichlet = dirichlet_values(model, X, Y)
UDIRICHLET = DIRICHLET_VALUES([X],[Y], MODEL) Examples dirichlet_values([0,1,2],[1,1,1],struct(name_dirichlet_values, homogeneous, ... c_dir, 1)) dirichlet_values([0:0.1:1],[0],struct(name_dirichlet_values, xstripes, ... c_dir, [0 1 2], ... dir_borders, [0.3 0.6])) dirichlet_values([0:0.1:1],[0],struct(name_dirichlet_values, box, ... c_dir, 1, ... dir_box_xrange, [0.3 0.6], ... dir_box_yrange, [-0.1 0.1]))
function [ L_E_neu , b_E_neu ] = fv_operators_neumann_explicit(model, model_data, U, NU_ind)
computes a neumann contribution matrix for finite volume time evolution operators, or their Frechet derivative
function rb_sim_data = rb_reconstruction_default(model, detailed_data, rb_sim_data)
(trivial) function computing a detailed reconstruction by linear combination of the coefficients in t...
function p = fv_plot(gridbase grid, dofs, params)
routine plotting a single fv function of fv_functions.
function a0 = rb_init_values_default(model, detailed_data)
function computing the reduced basis initial values. If the decomposition mode is coefficients...
function [ v , l2norm ] = fv_operators_output(model, model_data)
function returning components, coefficients, and complete operator for a linear output functional on ...
function p = plot_sequence(varargin)
plotting a sequence of data slices on polygonal 2d grid (constructed from params if empty) and provid...
function save_detailed_simulations(model, model_data, M, savepath)
perform loop over detailed simulations and save results or check consistency with existing saved resu...
function Umean = fv_element_mean(model, model_data, U, I)
function computing the element averages of a discrete function U in the grid elements with indices I...
function reduced_data_subset = lin_evol_reduced_data_subset(model, reduced_data)
method which modifies reduced_data, which is the data, that will be passed to the online-simulation a...
function [ flux , lambda ] = conv_flux_linear(glob, U, params)
function computing the convective flux of a convection problem.
function [ RBext , dummy ] = RB_extension_PCA_fixspace_flexible(model, detailed_data)
function computing a RB basis extension for given parameters by the POD-Greedy algorithm.