Skip to content
Snippets Groups Projects
Commit 20c88b92 authored by jcandrsn's avatar jcandrsn
Browse files

Important correction to adjoint indexing in uq_coll_add_response_adjoint. ...

Important correction to adjoint indexing in uq_coll_add_response_adjoint.  Corrections to boundary conditions for nonlinear vibration absorber vector field.  Update nonlinear robust optimization demo that uses the boundary condition corrections.
parent edeed07c
No related branches found
No related tags found
1 merge request!4Development
Showing
with 29 additions and 62 deletions
addpath('../../nonlinear_absorber_w_energy_harvester')
% addpath('../../utils')
% Initialize problem instance and set options
prob_init = coco_prob();
prob_init = coco_set(prob_init, 'coll', 'NTST', 15, 'NCOL', 4);
prob_init = coco_set(prob_init, 'cont', 'h', 0.5, 'h_max', 1000);
prob_init = coco_set(prob_init, 'cont', 'h', 0.5);
% prob_init = coco_set(prob_init, 'cont', 'PtMX', 100, 'NPR', 1, ...
% 'h', 0.5, 'almax', 30, 'h_max', 1000);
% prob_init = coco_set(prob_init, 'cont', 'NAdapt', 5, 'NPR', 10);
% % Temporarily allow really bad results for some checking of
% % Jacobians.
% prob_init = coco_set(prob_init, 'coll', 'TOL', 5);
run_name = 'uq_setup';
% Parameters
......@@ -68,7 +63,8 @@ if not(coco_exist(run_name, 'run'))
prob = coco_add_event(prob, 'VMAX', 'v0', 0);
bd = coco(prob, run_name, [], 1, {'v0', 'R', 'Om'});
% Phase shift until v0 = 0
bd = coco(prob, run_name, [], 1, {'v0', 'R'});
else
bd = coco_bd_read(run_name);
end
......@@ -80,7 +76,7 @@ if not(coco_exist(run_name, 'run')) || true
labs = coco_bd_labs(bd, 'VMAX');
prob = coco_set(prob_init, 'coll', 'NTST', 15, 'NCOL', 4);
prob = coco_set(prob, 'cont', 'PtMX', 300);
prob = coco_set(prob, 'cont', 'PtMX', 300, 'h_max', 1000);
prob = coco_set(prob, 'cont', 'NAdapt', 0, 'NPR', 10);
sol = coll_read_solution('bvp.seg1', 'uq_setup', labs(2));
......@@ -93,12 +89,11 @@ if not(coco_exist(run_name, 'run')) || true
% 'uq'-toolbox arguments
uq_args = {
{'z2', 'al'}, ... % Stochastic Parameter
{'Uniform', 'Uniform'}, ... % Distribution Type
[[0.01, 0.5]; [0.05, 0.5]], ... % [Mean, St.Dev]
{'Om'}, ...
'-add-adjt'}; % Denote Om as a unique
% deterministic parameter
{'z2'}, ... % Stochastic Parameter
{'Uniform'}, ... % Distribution Type
[[0.01, 0.5]], ... % [Mean, St.Dev]
{},...
'-add-adjt'};
prob = uq_isol2bvp_sample(prob, 'orig', ...
coll_args{:}, pnames, bvp_args{:}, uq_args{:});
......@@ -144,7 +139,7 @@ if not(coco_exist(run_name, 'run')) || true
{'rdo', 'b', 'd.rdo', 'orig.rmax.mean', 'orig.rmax.var', ...
'orig.vrms.mean', 'orig.vrms.var', 'd.up.z2', ...
'd.lo.z2', 'd.X0', 'd.z1', 'd.mu', ...
'd.Xi', 'd.la', 'd.ka', 'd.up.al', 'd.lo.al'}, [1.39, 10]);
'd.Xi', 'd.la', 'd.ka', 'd.al', 'd.Om'}, [1.39, 10]);
else
bd = coco_bd_read(run_name);
......
function r = x10(data, T0, T, x0, x1, p)
% Need to think through a helper function where the
% Jacobian is provided for a single sample and the helper
% function figures out how to blow it up into a Jacobian
% for all of the samples.
r = x0(1,:)./p(1,:);
end
\ No newline at end of file
......@@ -2,6 +2,6 @@ function y = nla_bc(data, T, x0, x1, p)
Om = p(10);
y = [x1(1:5)-x0(1:5); x1(6)-x0(6)-2*pi; T - 2*pi/Om];
y = [x1(1:5)-x0(1:5); x1(6)-x0(6)-2*pi];
end
......@@ -2,7 +2,7 @@ function [data, y] = nla_bc_caf(prob, data, u)
% Version that's compatible with coco_add_func
pr = data.pr;
bc = pr.bvp_bc;
bc = pr.bc_data;
x0 = u(bc.x0_idx);
x1 = u(bc.x1_idx);
......
......@@ -8,8 +8,6 @@ J = [ 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0;%-2*pi/Om^2];%;
1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -2*pi/Om^2;
];
0 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0;];%-2*pi/Om^2];%;
end
\ No newline at end of file
function [data, J] = nla_bc_du_caf(prob, data, u)
% T0 T x01 x02 x03 x04 x05 x06 x11 x12 x13 x14 x15 x16 X0 z1 mu z2 b al Xi la ka Om
J = [ 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0;
% 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
% T x01 x02 x03 x04 x05 x06 x11 x12 x13 x14 x15 x16 X0 z1 mu z2 b al Xi la ka Om
J = [ 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0;
];
end
\ No newline at end of file
function dJ = nla_bc_dudu(data, T, x0, x1, p)
Om = p(10);
% dJ = zeros(6,23,23);
dJ = zeros(7,23,23);
dJ(7,23,23) = 2*pi/Om^3;
dJ = zeros(6,23,23);
end
\ No newline at end of file
function y = nla_bc_v0(data, T, x0, x1, p)
Om = p(10);
y = [x1(1:5)-x0(1:5); x1(6)-x0(6)-2*pi; x0(2)];
y = [x1(1:5)-x0(1:5); x1(6)-x0(6)-2*pi; T - 2*pi/Om; x0(2)];
end
end
\ No newline at end of file
function J = nla_bc_v0_du(data, T, x0, x1, p)
Om = p(10);
% T x01 x02 x03 x04 x05 x06 x11 x12 x13 x14 x15 x16 X0 z1 mu z2 b al Xi la ka Om
J = [ 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
......@@ -9,8 +7,6 @@ J = [ 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 -1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0;
1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -2*pi/Om^2;
0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
];
0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;];
end
\ No newline at end of file
function dJ = nla_bc_v0_dudu(data, T, x0, x1, p)
Om = p(10);
dJ = zeros(8,23,23);
dJ(7,23,23) = 2*pi/Om^3;
dJ = zeros(7,23,23);
end
\ No newline at end of file
......@@ -62,9 +62,9 @@ elseif strcmp(uq_data.response_type, 'bv')
pend = 2+2*xdim+pdim;
for i=1:nsamples
[axidx, fdata] = coco_get_adjt_data(prob, ...
[axidx, adata] = coco_get_adjt_data(prob, ...
coco_get_id(uq_data.sids{i}, 'coll'), 'axidx', 'data');
maps = fdata.coll_seg.maps;
maps = adata.coll_opt;
aidx(i,1) = axidx(maps.T0_idx);
aidx(i,2) = axidx(maps.T_idx);
......@@ -100,11 +100,6 @@ if isfield(args, 'run')
mean_par_l0 = mean_par_chart.x;
mean_par_tl0 = cdata.v(mean_par_lidx);
else
% Will only work for a single response function. If
% adapting for higher dimension response functions,
% need to replace with zeros(num_response_outputs).
% num_response_outputs is not currently determined
% anywhere in the code.
mean_l0 = 0;
mean_tl0 = 0;
mean_par_l0 = 0;
......@@ -136,11 +131,6 @@ if isfield(args, 'run')
var_par_l0 = var_par_chart.x;
var_par_tl0 = cdata.v(var_par_lidx);
else
% Will only work for a single response function. If
% adapting for higher dimension response functions,
% need to replace with zeros(num_response_outputs).
% num_response_outputs is not currently determined
% anywhere in the code.
var_l0 = 0;
var_tl0 = 0;
var_par_l0 = 0;
......
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