1 function U0 = init_values_old_to_be_split_in_files(X,Y, params)
2 %
function U0 = init_values([X],[Y], params)
4 %
function constructing the initial values of the convection diffusion
5 % problem in the specified points and parameters
7 %
'blobs': u_0 = beta*exp(-gamma*((x-0.25).^2+(y-0.5).^2))* Xi(B_r) ...
8 % + (1-beta)*exp(-gamma*((x-0.25).^2+(y-0.7).^2))* Xi(B_r)
9 %
'homogeneous' u_0 = beta
10 %
'as_dirichlet' : evaluate dirichlet-data-
function at time 0
11 %
'waveproduct' : product of two sinus waves in both coordinate directions
12 %
'function_ptr' : init
function given by
function
13 % pointer to complete evaluation:
14 % params.init_values_ptr
16 % allow call U = f(X,Y,params), where X, Y are
17 % vectors of the same size. Return is a vector of
18 % values in the corresponding points.
19 %
'decomp_function_ptr' : init
function given by
function
20 % pointer to components and coefficients:
21 % params.init_values_coefficients_ptr
23 % allow call U = f([],[],params)
24 % Return is a vector of length Q with (parameter
25 % dependent) coefficients in U.
26 % params.init_values_components_ptr
28 % allow call U = f(X,Y,params), where X, Y are
29 % vectors of the same size.
30 % Return of components is a cell array of vectors of
31 % the same size as X and Y with the point values
32 % of the components. Linear combination of
33 % components by coefficients then yields the
34 % complete evaluation.
36 % required fields of params:
37 % name_init_values :
'blobs',
'homogeneous',
'wave',
'leftright',
40 %
if name_init_values ==
'gradient_box'
41 % c_init_up : constant value on upper border of box
43 % c_init_lo : constant value on lower border of box
45 %
if name_init_values ==
'blobs'
46 % radius : radius of blob-cutoff circle
47 % gamma : spread of initial-value gaussian
48 % beta : value between 0 and 1 weighting two gauss-blobs
50 %
if name_init_values ==
'homogeneous'
51 % c_init : constant value in field
53 %
if name_init_values ==
'wave'
54 % a wave of values between 0 and c_init
55 % c_init* 0.5 * (sin(freq_x * x + freq_y * y) + 1)
56 % c_init : maximum value in field
57 % freq_x : frequency in x-direction
58 % freq_y : frequency in y-direction
60 %
if name_init_values ==
'waveproduct'
61 % a product of axis-dependent waves of values between
62 % c_init_min and c_init_max
63 % c_init_max : maximum value in field
64 % c_init_min : minimum value in field
65 % c_init_freq_x : frequency in x-direction
66 % c_init_freq_y : frequency in y-direction
67 % c_init_phase_x : phase shift
68 % c_init_phase_y : phase shift
70 %
if name_init_values ==
'leftright'
71 % c_init_left : constant value in field left of border
72 % c_init_right : constant value in field left of border
73 % c_init_middle : x-coordinate separation
75 %
if name_init_values ==
'as_dirichlet'
76 % parameters as required by dirichlet
function
78 % Function supports affine decomposition, i.e. different operation modes
79 % guided by optional field decomp_mode in params. See also the
80 % contents.txt
for general explanation
82 % optional fields of params:
83 % mu_names : names of fields to be regarded as parameters in vector mu
85 % in
'coefficients' mode, the parameters in brackets are empty
87 % Bernard Haasdonk 5.10.2005
89 % determine affine_decomposition_mode as integer
90 decomp_mode = params.decomp_mode;
91 % flag indicating whether the computation respected the decomposition
92 respected_decomp_mode = 0;
94 if isequal(params.name_init_values,
'blobs')
95 B1 = ((X(:)-0.25).^2+(Y(:)-0.5).^2);
96 B2 = ((X(:)-0.25).^2+(Y(:)-0.7).^2);
97 I1 = (B1<params.radius.^2);
98 I2 = (B2<params.radius.^2);
100 U0 = params.beta*exp(-params.gamma*B1).*I1 ...
101 + (1-params.beta)* exp(-params.gamma*B2).*I2;
102 elseif decomp_mode == 1
103 if ismember(
'gamma',params.mu_names) | ...
104 ismember(
'radius',params.mu_names)
105 error('affine decomp with respect to mu_names not possible!');
108 U0{1} = exp(-params.gamma*B1).*I1;
109 U0{2} = exp(-params.gamma*B2).*I2;
110 else % decomp_mode == 2
111 U0 = [params.beta, 1-params.beta];
113 respected_decomp_mode = 1;
115 elseif isequal(params.name_init_values,
'leftright')
116 i = (X <= params.init_middle);
117 U0 = params.c_init_left * i + (1-i) * params.c_init_right;
119 elseif isequal(params.name_init_values,'decomp_function_ptr')
121 U0 = params.init_values_coefficients_ptr([],[],params);
122 elseif decomp_mode == 1; % components
123 U0 = params.init_values_components_ptr(X,Y,params);
124 else % decomp_mode = 0, complete
125 Ucoefficients = params.init_values_coefficients_ptr([],[],params);
126 Ucomponents = params.init_values_components_ptr(X,Y,params);
127 U0 = lincomb_sequence(Ucomponents,Ucoefficients);
129 respected_decomp_mode = 1;
131 elseif isequal(params.name_init_values,'function_ptr')
132 U0 = params.init_values_ptr(X,Y,params);
133 respected_decomp_mode = 0;
136 error('unknown name_init_values');
139 if decomp_mode>0 && respected_decomp_mode==0
140 error('function does not support affine decomposition!');
143 % TO BE ADJUSTED TO NEW SYNTAX