You can not select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
188 lines
6.8 KiB
188 lines
6.8 KiB
3 years ago
|
%--------------------------
|
||
|
% @Author: Jingqiao Hu
|
||
|
% @Date: 2021-08-25 17:26:51
|
||
|
% @LastEditTime: 2021-09-09 14:25:42
|
||
|
%--------------------------
|
||
|
function [f0val, df0dx, fval, dfdx] = obj_cons_sens_micro(nelx, nely, microx, H_global, Hs_global, fext, U, nele_global, v0_micro, ...
|
||
|
dk_micro, edofMat_ma, rho_cell, cell2global_id, e0, e1, penal, ke1, iK_mi, jK_mi, edofMat_mi, bdofs_num, dD, dofid, newdofs, ...
|
||
|
bdofs_plus, ib_dofs, ii_dofs, regulated_matrix, bezier_B, dFdx, De)
|
||
|
|
||
|
pnorm = 16;
|
||
|
nele_macro = nelx * nely;
|
||
|
nele_micro = microx^2;
|
||
|
|
||
|
%% compliance & volume_macro
|
||
|
comp = fext' * U;
|
||
|
|
||
|
% compliance_sensitivities(nelx, nely, microx, dk_micro, U, edofMat_ma, 1, 1);
|
||
|
|
||
|
% pnorm-volume constraint of micro
|
||
|
volume_micro = reshape(mean(rho_cell, 1), [], 1); % [m,M] -> [1,M]
|
||
|
cons_v_micro = (mean(volume_micro.^pnorm))^(1/pnorm) / v0_micro - 1;
|
||
|
|
||
|
dpnorm = 1/nele_macro * (mean(volume_micro.^pnorm))^(1/pnorm-1) * (volume_micro.^(pnorm-1)) / v0_micro;
|
||
|
dpn_expand = kron(ones(nele_micro, 1), dpnorm); % expand to [nele_global, 1]
|
||
|
dv_micro = ones(nele_global, 1) / nele_global .* dpn_expand;
|
||
|
|
||
|
% cons_v_micro = mean(rho_cell(:)) / v0_micro - 1;
|
||
|
% dv_micro = ones(nele_global, 1) / nele_global / v0_micro;
|
||
|
|
||
|
[dc_micro] = obj_sens_micro(nelx, nely, microx, dk_micro, U, edofMat_ma, cell2global_id);
|
||
|
|
||
|
dc_micro1 = compliance_sensitivities(nelx, nely, microx, regulated_matrix, U, edofMat_ma, edofMat_mi, dofid, newdofs, ...
|
||
|
e0, e1, penal, ke1, iK_mi, jK_mi, bdofs_num, dD, rho_cell, bdofs_plus, ib_dofs, ii_dofs, cell2global_id, bezier_B, dFdx, De);
|
||
|
|
||
|
max(dc_micro(:) - dc_micro1(:))
|
||
|
|
||
|
dc_micro(:) = H_global * (dc_micro(:) ./ Hs_global);
|
||
|
dv_micro(:) = H_global * (dv_micro(:) ./ Hs_global );
|
||
|
% CHECK: maybe need to be scaled to 1000 times
|
||
|
|
||
|
%% sensitivities
|
||
|
f0val = comp;
|
||
|
df0dx = - dc_micro(:); % note the negetive sign
|
||
|
|
||
|
fval = cons_v_micro;
|
||
|
dfdx = dv_micro(:)';
|
||
|
end
|
||
|
|
||
|
% newdofs = [boundarys, inner];
|
||
|
function dc_micro = compliance_sensitivities(nelx, nely, microx, regulated_matrix, U, edofMat_ma, edofMat_mi, dofid, newdofs, ...
|
||
|
e0, e1, penal, ke1, iK_mi, jK_mi, bdofs_num, dD, rho, bdofs_plus, ib_dofs, ii_dofs, cell2global_id, bezier_B, dFdx, De)
|
||
|
|
||
|
nele_macro = nelx * nely;
|
||
|
nele_micro = microx^2;
|
||
|
UM = U(edofMat_ma); % [M,v]
|
||
|
|
||
|
order1 = dofid(edofMat_mi'); % order as [border,idofs]
|
||
|
order2 = node2ele_id(dofid(edofMat_mi));
|
||
|
|
||
|
% O2 = zeros(2*(microx+1)^2 - bdofs_num, 1);
|
||
|
O = zeros(bdofs_num, 1);
|
||
|
idofs_plus = 1 - bdofs_plus; % [8,m]
|
||
|
|
||
|
dc = zeros(nele_micro, nele_macro);
|
||
|
dc1 = zeros(nele_micro, nele_macro);
|
||
|
|
||
|
for ele = 1: nele_macro
|
||
|
rho_e = rho(:, ele);
|
||
|
rho_penal = e0 + rho_e.^penal * (e1-e0); % m,1 micro-rho need to be penal
|
||
|
drho = penal * rho_e.^(penal-1) * (e1 - e0);
|
||
|
|
||
|
ue = UM(ele, :); % [1,v]
|
||
|
Re = regulated_matrix{ele}; % [2n,v]
|
||
|
|
||
|
former1 = Re * ue'; % [2n,1]
|
||
|
former3 = former1(order1); % expand to [8,m]
|
||
|
|
||
|
former4 = former3' * dD .* rho_penal; % [m,8]
|
||
|
former5 = reshape(former4', [], 1);
|
||
|
former6 = former5(order2); % assemble to [2n,1]
|
||
|
former_sub = former6(bdofs_num + 1 : end); % choose dofs corresponding to k_22 of [2i,1]
|
||
|
|
||
|
K = assemble_micro_k(dofid, rho_e, microx, e0, e1, penal, ke1, iK_mi, jK_mi);
|
||
|
k22 = K(bdofs_num + 1 : end, bdofs_num + 1 : end); % [2i,2i]
|
||
|
|
||
|
former_sub1 = k22 \ former_sub; % [2i,1]
|
||
|
former7 = [O; former_sub1]'; % [1,2n]
|
||
|
former = former7(order1); % [8,m]
|
||
|
|
||
|
latter_b(:,1,:) = former3 .* bdofs_plus; % [8,1,m]
|
||
|
latter_i(:,1,:) = - former3 .* idofs_plus; % [8,1,m]
|
||
|
|
||
|
dke = ke1(:) * drho(:)'; % [8*8, m]
|
||
|
dke1 = reshape(dke, 8, 8, nele_micro); % [8, 8, m]
|
||
|
|
||
|
dk21 = dke1 .* ib_dofs; % [8, 8, m]
|
||
|
dk22 = dke1 .* ii_dofs; % [8, 8, m]
|
||
|
|
||
|
dM = mtimesx(dk22, latter_i) - mtimesx(dk21, latter_b); % [8,1,m]
|
||
|
dM1 = squeeze(dM); % [8,m]
|
||
|
|
||
|
% dM1 = reshape(dM1, [], 1);
|
||
|
% dM2 = dM1(order2); % assemble to [2n,1]
|
||
|
% dM2_sub = dM2(bdofs_num + 1 : end); % [2i,1]
|
||
|
% latter_sub = k22 \ dM2_sub; % [2i,1]
|
||
|
%
|
||
|
% latter = [O; latter_sub];
|
||
|
% latter = latter(order1); % [8,m]
|
||
|
|
||
|
dc(:,ele) = 2 * sum(former .* dM1, 1); % [1,m]
|
||
|
% dc(:,ele) = 2 * sum(former4' .* latter, 1); % [1,m]
|
||
|
|
||
|
%%
|
||
|
ndofs = 2 * (microx + 1)^2;
|
||
|
dk = zeros(ndofs, ndofs);
|
||
|
O = zeros(bdofs_num, 1);
|
||
|
for ele_mi = 12:nele_micro
|
||
|
edof = dofid(edofMat_mi(ele_mi, :));
|
||
|
dke = ke1 * drho(ele_mi);
|
||
|
dk(edof, edof) = dke; % assemble by the order of [dofid]
|
||
|
|
||
|
dk21 = dk(bdofs_num+1 : end, 1 : bdofs_num) * former1(1:bdofs_num);
|
||
|
dk22 = dk(bdofs_num+1 : end, bdofs_num+1 : end) * (-former1(bdofs_num+1:end));
|
||
|
dk2_pre = k22 * (dk22 - dk21); % NOTE: not including \bO
|
||
|
|
||
|
dphi = [O; dk2_pre];
|
||
|
% a = dphi(edof);
|
||
|
%
|
||
|
% test = ke1 \ dM1(:,ele_mi);
|
||
|
|
||
|
% a = former1(edof)' * dD * rho_penal(ele_mi);
|
||
|
% max(dM1(:,ele_mi) - a)
|
||
|
|
||
|
dc1(ele_mi,ele) = 2 * former1(edof)' * dD * rho_penal(ele_mi) * dphi(edof); % + former1(edof)' * dD * drho(ele_mi) * former1(edof);
|
||
|
% dc2(ele_mi,ele) = 2 * former6(edof)' * dphi(edof);
|
||
|
dk(:) = 0;
|
||
|
end
|
||
|
end
|
||
|
|
||
|
dc_micro = dc(cell2global_id);
|
||
|
dc_micro1 = dc1(cell2global_id);
|
||
|
end
|
||
|
|
||
|
function order_id = node2ele_id(edofMat_mi)
|
||
|
order1 = reshape(edofMat_mi', [], 1);
|
||
|
order1(:,2) = 1:size(order1,1);
|
||
|
order2 = sortrows(order1, 1);
|
||
|
[~,ia,ic] = unique(order2(:,1)); % C = A(ia)
|
||
|
idx = order2(:,2);
|
||
|
order_id = idx(ia);
|
||
|
end
|
||
|
|
||
|
% compliance_micro sensitivities
|
||
|
function [dc_micro] = obj_sens_micro(nelx, nely, microx, dk_micro, U, edofMat_ma, cell2global_id)
|
||
|
|
||
|
nele_micro = microx^2;
|
||
|
nele_macro = nelx * nely;
|
||
|
vdofs_num = size(edofMat_ma, 2);
|
||
|
|
||
|
% tic
|
||
|
UM = U(edofMat_ma); % [M,v]
|
||
|
UM1 = kron(ones(nele_micro,1), UM)'; % [v,m*M]
|
||
|
UM2 = reshape(UM1, vdofs_num, 1, nele_macro, nele_micro); % [v,1,M,m]
|
||
|
UM3 = permute(UM2, [1,2,4,3]); % [v,1,m,M]
|
||
|
UM4 = permute(UM3, [2,1,3,4]); % [1,v,m,M]
|
||
|
|
||
|
dc5 = squeeze(mtimesx(mtimesx(UM4, dk_micro), UM3)); % [m,M]
|
||
|
dc_micro = dc5(cell2global_id);
|
||
|
% toc
|
||
|
|
||
|
% tic
|
||
|
% dc = zeros(nele_micro, nele_macro);
|
||
|
% parfor ele = 1: nelx * nely
|
||
|
% ue = UM(ele, :);
|
||
|
%
|
||
|
% dke = dk_micro(:,:,:,ele); % [v,v,m]
|
||
|
% dk1 = reshape(dke, vdofs_num, []); % [v,v*m]
|
||
|
%
|
||
|
% dc_former = ue * dk1; % % [1,v*m]
|
||
|
% dc1 = reshape(dc_former, vdofs_num, nele_micro);
|
||
|
%
|
||
|
% dc(:,ele) = dc1' * ue';
|
||
|
% end
|
||
|
% dc_micro1 = dc(cell2global_id);
|
||
|
% toc
|
||
|
%
|
||
|
% max(dc_micro(:) - dc_micro1(:))
|
||
|
end
|