이 질문을 팔로우합니다.
- 팔로우하는 게시물 피드에서 업데이트를 확인할 수 있습니다.
- 정보 수신 기본 설정에 따라 이메일을 받을 수 있습니다.
I wrote code on MATLAB to solve Stokes equations using Bramble and Pasciak method but I got strange values of error for pressure p
조회 수: 3 (최근 30일)
이전 댓글 표시
%Define problem parameters
mu = 1.0; % viscosity
%f = [1; 0; 1; 2]; % forcing term for velocity
%g = [2; 1]; % forcing term for pressure
%% Mesh creation
h = 0.5; % discretization step
[node,elem] = squaremesh([-1,1,-1,1],h);
Unrecognized function or variable 'squaremesh'.
subplot(1,2,1); showmesh(node,elem);
N = size(node,1);
NTc = size(elem,1);
% Uniform refinement to get a fine grid
[node,elem] = uniformrefine(node,elem);
% Uniform refinement to get a fine grid
subplot(1,2,2); showmesh(node,elem);
%%extract boundary conditions
[Dirichlet, bdFlag] = ExtractBoundaryEdges(elem, node);
%% Load Stokes Data Partial Differential Equations
%Different velocity and pressure equation models
pde = StokesModelMustafaNew;%StokesModelMustafa;%StokesModelMustafaNew;%Stokesdata3;%StokesModelMustafa2;%StokesModelMustafaSin;%;
option.elemType = 'isoP2P0';
option.fquadorder = 3;
tolerance = 1e-3;
max_iterations = 1000;
mesh.node = node;
mesh.elem = elem;
mesh.bdFlag = bdFlag;
%The following line of code is only run to get the eqn (equation structure) easier without defining it every time
[soln,eqn,info] = StokesisoP2P0(node,elem,bdFlag,pde, option);
%% Call the solvers and calulate L2 errors, run time and outcome as well as convergence rate
% [soln,eqn,err,time,solver] = femStokes(mesh,pde,option);
% return
uI = pde.exactu([node; (node(eqn.edge(:,1),:)+node(eqn.edge(:,2),:))/2]);
uI = reshape(uI, [2*length(uI),1]);
pI = pde.exactp([node(elem(:,1),1) node(elem(:,3),2)]);
% Define matrix H
%H = [1 0 0 2; 0 2 1 0; 2 3 1 0; 0 1 0 1];
% Define matrix A
%A = H * H';
% Define matrix B
%B = [1 2 1 0; 0 1 0 0];
% Assemble the finite element system
%[A, B, ~, ~] = assempde(pde,node,elem,'f',f,'g',g,'mu',mu);
tic;
% Assemble the Stokes matrix
N = size(eqn.A,1);
M = size(eqn.B,1);
%StokesMatrix = [eqn.A, eqn.B'; eqn.B, sparse(M,M)];
% Define dimensions of eqn.B'
eqn.B_transpose = eqn.B'; % Transpose of eqn.B
% Construct a sparse MxM matrix (adjust M to your desired size)
M = size(eqn.B, 1); % Number of columns in eqn.B (or rows in eqn.B')
sparse_M_M = sparse(M, M);
% Check dimensions of matrices involved in the problematic line
size_eqn.B = size(eqn.B);
size_eqn.A = size(eqn.A);
% Perform the matrix multiplication explicitly and verify dimensions
result = 2 * eqn.B / (eqn.A)* eqn.B_transpose;
size_result = size(result);
% Concatenate eqn.A, eqn.B_transpose, eqn.B, and sparse_M_M
StokesMatrix = [eqn.A, eqn.B_transpose; eqn.B, sparse_M_M];
% Define matrix M
%M = [2*eye(size(eqn.A,1)) 2*(eqn.A)\eqn.B_transpose; eqn.B 2*eqn.B*(eqn.A)\eqn.B_transpose];
M_top = [2 * eye(size(eqn.A, 1)), 2 * (eqn.A\eqn.B_transpose)];
M_bottom = [eqn.B, 2 * eqn.B * (eqn.A \ eqn.B_transpose)];
M = [M_top; M_bottom];
% Define matrix F
F = [2 * eqn.A \ eqn.f; 2 * eqn.B * (eqn.A \ eqn.f) - eqn.g];
% Define the right-hand side
rhs = [eqn.f; eqn.g];
x0 = ones(size(F));
% Solve the Stokes system
solution = StokesMatrix \ rhs;
sol = M \ F;
% Extract the velocity and pressure components from the solution
u = solution(1:N);
p = solution(N+1:end);
% Calculate the residual R
R = F - M * sol;
% Define matrix C
C = [0.5 * eqn.A zeros(size(transpose(eqn.B))); zeros(size(eqn.B)) eye(size(eqn.B,1))];
% Define an appropriate non-zero vector x
[res_conjgrad, iterCg, relresCg, H1ErrUCg, L2ErrPCg, L2ErrGradUandPCg] = conjgradNew3(eqn.A, C,F ,[eqn.f; eqn.g],x0 , max_iterations, tolerance, u, p, h, M);
for i = 1:iterCg
if H1ErrUCg(i) == 0 && i == 1
H1ErrUCg(i) = max(H1ErrUCg);
elseif H1ErrUCg(i) == 0
H1ErrUCg(i) = H1ErrUCg(i-1);
end
end
for i = 1:iterCg
if L2ErrPCg(i) == 0 && i == 1
L2ErrPCg(i) = max(L2ErrPCg);
elseif L2ErrPCg(i) == 0
L2ErrPCg(i) = L2ErrPCg(i-1);
end
end
for i = 1:iterCg
if L2ErrGradUandPCg(i) == 0 && i == 1
L2ErrGradUandPCg(i) = max(L2ErrGradUandPCg);
elseif L2ErrGradUandPCg(i) == 0
L2ErrGradUandPCg(i) = L2ErrGradUandPCg(i-1);
end
end
p = res_conjgrad(size(eqn.A,1)+1:end);
u = res_conjgrad(1:size(eqn.A,1));
time_BPCG_finish = toc;
function [inner_pro] = inner_x_y(C, x,y)
%C = [A 0 ; 0 1]
%(x,y) = (Cx,y) = (Ax_1,y_1) + (x_2,y_2) = transpose(x_1)*A*y_1 +transpose(x_2)*y_2
%trans(y_1)*x_2
inner_pro = transpose(y) * C * x;
end
댓글 수: 17
Mostafa Kassoum
2024년 1월 12일
편집: Walter Roberson
2024년 1월 12일
The sequare mesh comes from the following link;
Mostafa Kassoum
2024년 1월 12일
편집: Walter Roberson
2024년 1월 12일
Also to test the code you need the following code:
function [x, i, relres, H1ErrU, L2ErrP, L2ErrGradUandP] = conjgradNew3(~, C,F,~, x, maxiter, tol, u, p, h, M)
z_old = x;
p_old = F - M*z_old;
r_old=p_old;
for i = 1:maxiter
alpha_i = inner_x_y(C, r_old,p_old) / inner_x_y(C, M*p_old,p_old);
z_new = z_old + alpha_i*p_old;
r_new = F - M*z_new;
norm_r = transpose(r_new)*r_new;%inner_x_y(C, r_new,r_new);
H1ErrU(i) = norm(gradient(u - z_new(1:length(u))));
L2ErrP(i) = norm(p - z_new(length(u)+1:end));
L2ErrGradUandP(i) = sqrt((H1ErrU(i))^2 + (L2ErrP(i))^2);
%L2ErrGradUandP(i) = sqrt((norm(gradient(z_new(1:length(uI)),h) - gradient(uI,h)))^2 + (norm(pI - z_new(length(uI)+1:end)))^2);
relres(i) = norm_r;
fprintf('#dof: %8.0u, Bramble Pasciak CG iter: %2.0u, err = %12.8g\n \n',...
length(z_new)+length(z_new(length(u)+1:end)), i, norm_r);
if (norm_r)< tol
break
end
beta_i = inner_x_y(C, M*r_new,p_old) / inner_x_y(C, M*p_old,p_old);
p_new = r_new - beta_i*p_old;
% Update
z_old = z_new;
r_old = r_new;
p_old = p_new;
end
end
Walter Roberson
2024년 1월 12일
The given archive does not include ExtractBoundaryEdges or StokesModelMustafaNew
%bring in external tools
sqm_url = "https://raw.githubusercontent.com/lyc102/ifem/master/mesh/squaremesh.m";
shm_url = "https://raw.githubusercontent.com/lyc102/ifem/master/tool/showmesh.m";
unir_url = "https://raw.githubusercontent.com/lyc102/ifem/master/mesh/uniformrefine.m";
myu_url = "https://raw.githubusercontent.com/lyc102/ifem/master/tool/myunique.m";
tn = tempname();
sqm_localname = fullfile(tn, "squaremesh.m");
shm_localname = fullfile(tn, "showmesh.m");
unir_localname = fullfile(tn, "uniformrefine.m");
myu_localname = fullfile(tn, "myunique.m");
mkdir(tn);
addpath(tn)
websave(sqm_localname, sqm_url);
websave(shm_localname, shm_url);
websave(unir_localname, unir_url);
websave(myu_localname, myu_url);
%done bringing in external tools
%Define problem parameters
mu = 1.0; % viscosity
%f = [1; 0; 1; 2]; % forcing term for velocity
%g = [2; 1]; % forcing term for pressure
%% Mesh creation
h = 0.5; % discretization step
[node,elem] = squaremesh([-1,1,-1,1],h);
subplot(1,2,1); showmesh(node,elem);
N = size(node,1);
NTc = size(elem,1);
% Uniform refinement to get a fine grid
[node,elem] = uniformrefine(node,elem);
% Uniform refinement to get a fine grid
subplot(1,2,2); showmesh(node,elem);
%%extract boundary conditions
[Dirichlet, bdFlag] = ExtractBoundaryEdges(elem, node);
Unrecognized function or variable 'ExtractBoundaryEdges'.
%% Load Stokes Data Partial Differential Equations
%Different velocity and pressure equation models
pde = StokesModelMustafaNew;%StokesModelMustafa;%StokesModelMustafaNew;%Stokesdata3;%StokesModelMustafa2;%StokesModelMustafaSin;%;
option.elemType = 'isoP2P0';
option.fquadorder = 3;
tolerance = 1e-3;
max_iterations = 1000;
mesh.node = node;
mesh.elem = elem;
mesh.bdFlag = bdFlag;
%The following line of code is only run to get the eqn (equation structure) easier without defining it every time
[soln,eqn,info] = StokesisoP2P0(node,elem,bdFlag,pde, option);
%% Call the solvers and calulate L2 errors, run time and outcome as well as convergence rate
% [soln,eqn,err,time,solver] = femStokes(mesh,pde,option);
% return
uI = pde.exactu([node; (node(eqn.edge(:,1),:)+node(eqn.edge(:,2),:))/2]);
uI = reshape(uI, [2*length(uI),1]);
pI = pde.exactp([node(elem(:,1),1) node(elem(:,3),2)]);
% Define matrix H
%H = [1 0 0 2; 0 2 1 0; 2 3 1 0; 0 1 0 1];
% Define matrix A
%A = H * H';
% Define matrix B
%B = [1 2 1 0; 0 1 0 0];
% Assemble the finite element system
%[A, B, ~, ~] = assempde(pde,node,elem,'f',f,'g',g,'mu',mu);
tic;
% Assemble the Stokes matrix
N = size(eqn.A,1);
M = size(eqn.B,1);
%StokesMatrix = [eqn.A, eqn.B'; eqn.B, sparse(M,M)];
% Define dimensions of eqn.B'
eqn.B_transpose = eqn.B'; % Transpose of eqn.B
% Construct a sparse MxM matrix (adjust M to your desired size)
M = size(eqn.B, 1); % Number of columns in eqn.B (or rows in eqn.B')
sparse_M_M = sparse(M, M);
% Check dimensions of matrices involved in the problematic line
size_eqn.B = size(eqn.B);
size_eqn.A = size(eqn.A);
% Perform the matrix multiplication explicitly and verify dimensions
result = 2 * eqn.B / (eqn.A)* eqn.B_transpose;
size_result = size(result);
% Concatenate eqn.A, eqn.B_transpose, eqn.B, and sparse_M_M
StokesMatrix = [eqn.A, eqn.B_transpose; eqn.B, sparse_M_M];
% Define matrix M
%M = [2*eye(size(eqn.A,1)) 2*(eqn.A)\eqn.B_transpose; eqn.B 2*eqn.B*(eqn.A)\eqn.B_transpose];
M_top = [2 * eye(size(eqn.A, 1)), 2 * (eqn.A\eqn.B_transpose)];
M_bottom = [eqn.B, 2 * eqn.B * (eqn.A \ eqn.B_transpose)];
M = [M_top; M_bottom];
% Define matrix F
F = [2 * eqn.A \ eqn.f; 2 * eqn.B * (eqn.A \ eqn.f) - eqn.g];
% Define the right-hand side
rhs = [eqn.f; eqn.g];
x0 = ones(size(F));
% Solve the Stokes system
solution = StokesMatrix \ rhs;
sol = M \ F;
% Extract the velocity and pressure components from the solution
u = solution(1:N);
p = solution(N+1:end);
% Calculate the residual R
R = F - M * sol;
% Define matrix C
C = [0.5 * eqn.A zeros(size(transpose(eqn.B))); zeros(size(eqn.B)) eye(size(eqn.B,1))];
% Define an appropriate non-zero vector x
[res_conjgrad, iterCg, relresCg, H1ErrUCg, L2ErrPCg, L2ErrGradUandPCg] = conjgradNew3(eqn.A, C,F ,[eqn.f; eqn.g],x0 , max_iterations, tolerance, u, p, h, M);
for i = 1:iterCg
if H1ErrUCg(i) == 0 && i == 1
H1ErrUCg(i) = max(H1ErrUCg);
elseif H1ErrUCg(i) == 0
H1ErrUCg(i) = H1ErrUCg(i-1);
end
end
for i = 1:iterCg
if L2ErrPCg(i) == 0 && i == 1
L2ErrPCg(i) = max(L2ErrPCg);
elseif L2ErrPCg(i) == 0
L2ErrPCg(i) = L2ErrPCg(i-1);
end
end
for i = 1:iterCg
if L2ErrGradUandPCg(i) == 0 && i == 1
L2ErrGradUandPCg(i) = max(L2ErrGradUandPCg);
elseif L2ErrGradUandPCg(i) == 0
L2ErrGradUandPCg(i) = L2ErrGradUandPCg(i-1);
end
end
p = res_conjgrad(size(eqn.A,1)+1:end);
u = res_conjgrad(1:size(eqn.A,1));
time_BPCG_finish = toc;
function [inner_pro] = inner_x_y(C, x,y)
%C = [A 0 ; 0 1]
%(x,y) = (Cx,y) = (Ax_1,y_1) + (x_2,y_2) = transpose(x_1)*A*y_1 +transpose(x_2)*y_2
%trans(y_1)*x_2
inner_pro = transpose(y) * C * x;
end
function [x, i, relres, H1ErrU, L2ErrP, L2ErrGradUandP] = conjgradNew3(~, C,F,~, x, maxiter, tol, u, p, h, M)
z_old = x;
p_old = F - M*z_old;
r_old=p_old;
for i = 1:maxiter
alpha_i = inner_x_y(C, r_old,p_old) / inner_x_y(C, M*p_old,p_old);
z_new = z_old + alpha_i*p_old;
r_new = F - M*z_new;
norm_r = transpose(r_new)*r_new;%inner_x_y(C, r_new,r_new);
H1ErrU(i) = norm(gradient(u - z_new(1:length(u))));
L2ErrP(i) = norm(p - z_new(length(u)+1:end));
L2ErrGradUandP(i) = sqrt((H1ErrU(i))^2 + (L2ErrP(i))^2);
%L2ErrGradUandP(i) = sqrt((norm(gradient(z_new(1:length(uI)),h) - gradient(uI,h)))^2 + (norm(pI - z_new(length(uI)+1:end)))^2);
relres(i) = norm_r;
fprintf('#dof: %8.0u, Bramble Pasciak CG iter: %2.0u, err = %12.8g\n \n',...
length(z_new)+length(z_new(length(u)+1:end)), i, norm_r);
if (norm_r)< tol
break
end
beta_i = inner_x_y(C, M*r_new,p_old) / inner_x_y(C, M*p_old,p_old);
p_new = r_new - beta_i*p_old;
% Update
z_old = z_new;
r_old = r_new;
p_old = p_new;
end
end
Mostafa Kassoum
2024년 1월 12일
편집: Walter Roberson
2024년 1월 12일
The following code is the StokesModelMustafaNew
function pde = StokesModelMustafaNew
%% Custom Stokes equation model
%
pde = struct('f', 0, 'exactp', @exactp, 'exactu', @exactu,'g_D',@g_D,'g_N',@g_N);
% exact velocity
function z = exactu(p)
x = p(:,1);
y = p(:,2);
z(:,1) = 1/(5*pi()^2).*sin(pi()*x).*sin(2*pi()*y);
z(:,2) = 1/(5*pi()^2).*sin(2*pi()*x).*sin(pi()*y);
end
% exact pressure
function z = exactp(p)
x = p(:,1); y = p(:,2);
z = 2/3-x.^2-y.^2;
end
% Dirichlet boundary condition of velocity
function z = g_D(p)
z = exactu(p);
end
function z = g_N(p)
z(:,1) = 2*ones(size(p,1),1);
z(:,2) = zeros(size(p,1),1);
end
end
Mostafa Kassoum
2024년 1월 12일
I think now you can run the code and see why is the L^2 error of pressure p big value.
Walter Roberson
2024년 1월 12일
Still need ExtractBoundaryEdges
%bring in external tools
sqm_url = "https://raw.githubusercontent.com/lyc102/ifem/master/mesh/squaremesh.m";
shm_url = "https://raw.githubusercontent.com/lyc102/ifem/master/tool/showmesh.m";
unir_url = "https://raw.githubusercontent.com/lyc102/ifem/master/mesh/uniformrefine.m";
myu_url = "https://raw.githubusercontent.com/lyc102/ifem/master/tool/myunique.m";
tn = tempname();
sqm_localname = fullfile(tn, "squaremesh.m");
shm_localname = fullfile(tn, "showmesh.m");
unir_localname = fullfile(tn, "uniformrefine.m");
myu_localname = fullfile(tn, "myunique.m");
mkdir(tn);
addpath(tn)
websave(sqm_localname, sqm_url);
websave(shm_localname, shm_url);
websave(unir_localname, unir_url);
websave(myu_localname, myu_url);
%done bringing in external tools
%Define problem parameters
mu = 1.0; % viscosity
%f = [1; 0; 1; 2]; % forcing term for velocity
%g = [2; 1]; % forcing term for pressure
%% Mesh creation
h = 0.5; % discretization step
[node,elem] = squaremesh([-1,1,-1,1],h);
subplot(1,2,1); showmesh(node,elem);
N = size(node,1);
NTc = size(elem,1);
% Uniform refinement to get a fine grid
[node,elem] = uniformrefine(node,elem);
% Uniform refinement to get a fine grid
subplot(1,2,2); showmesh(node,elem);
%%extract boundary conditions
[Dirichlet, bdFlag] = ExtractBoundaryEdges(elem, node);
Unrecognized function or variable 'ExtractBoundaryEdges'.
%% Load Stokes Data Partial Differential Equations
%Different velocity and pressure equation models
pde = StokesModelMustafaNew;%StokesModelMustafa;%StokesModelMustafaNew;%Stokesdata3;%StokesModelMustafa2;%StokesModelMustafaSin;%;
option.elemType = 'isoP2P0';
option.fquadorder = 3;
tolerance = 1e-3;
max_iterations = 1000;
mesh.node = node;
mesh.elem = elem;
mesh.bdFlag = bdFlag;
%The following line of code is only run to get the eqn (equation structure) easier without defining it every time
[soln,eqn,info] = StokesisoP2P0(node,elem,bdFlag,pde, option);
%% Call the solvers and calulate L2 errors, run time and outcome as well as convergence rate
% [soln,eqn,err,time,solver] = femStokes(mesh,pde,option);
% return
uI = pde.exactu([node; (node(eqn.edge(:,1),:)+node(eqn.edge(:,2),:))/2]);
uI = reshape(uI, [2*length(uI),1]);
pI = pde.exactp([node(elem(:,1),1) node(elem(:,3),2)]);
% Define matrix H
%H = [1 0 0 2; 0 2 1 0; 2 3 1 0; 0 1 0 1];
% Define matrix A
%A = H * H';
% Define matrix B
%B = [1 2 1 0; 0 1 0 0];
% Assemble the finite element system
%[A, B, ~, ~] = assempde(pde,node,elem,'f',f,'g',g,'mu',mu);
tic;
% Assemble the Stokes matrix
N = size(eqn.A,1);
M = size(eqn.B,1);
%StokesMatrix = [eqn.A, eqn.B'; eqn.B, sparse(M,M)];
% Define dimensions of eqn.B'
eqn.B_transpose = eqn.B'; % Transpose of eqn.B
% Construct a sparse MxM matrix (adjust M to your desired size)
M = size(eqn.B, 1); % Number of columns in eqn.B (or rows in eqn.B')
sparse_M_M = sparse(M, M);
% Check dimensions of matrices involved in the problematic line
size_eqn.B = size(eqn.B);
size_eqn.A = size(eqn.A);
% Perform the matrix multiplication explicitly and verify dimensions
result = 2 * eqn.B / (eqn.A)* eqn.B_transpose;
size_result = size(result);
% Concatenate eqn.A, eqn.B_transpose, eqn.B, and sparse_M_M
StokesMatrix = [eqn.A, eqn.B_transpose; eqn.B, sparse_M_M];
% Define matrix M
%M = [2*eye(size(eqn.A,1)) 2*(eqn.A)\eqn.B_transpose; eqn.B 2*eqn.B*(eqn.A)\eqn.B_transpose];
M_top = [2 * eye(size(eqn.A, 1)), 2 * (eqn.A\eqn.B_transpose)];
M_bottom = [eqn.B, 2 * eqn.B * (eqn.A \ eqn.B_transpose)];
M = [M_top; M_bottom];
% Define matrix F
F = [2 * eqn.A \ eqn.f; 2 * eqn.B * (eqn.A \ eqn.f) - eqn.g];
% Define the right-hand side
rhs = [eqn.f; eqn.g];
x0 = ones(size(F));
% Solve the Stokes system
solution = StokesMatrix \ rhs;
sol = M \ F;
% Extract the velocity and pressure components from the solution
u = solution(1:N);
p = solution(N+1:end);
% Calculate the residual R
R = F - M * sol;
% Define matrix C
C = [0.5 * eqn.A zeros(size(transpose(eqn.B))); zeros(size(eqn.B)) eye(size(eqn.B,1))];
% Define an appropriate non-zero vector x
[res_conjgrad, iterCg, relresCg, H1ErrUCg, L2ErrPCg, L2ErrGradUandPCg] = conjgradNew3(eqn.A, C,F ,[eqn.f; eqn.g],x0 , max_iterations, tolerance, u, p, h, M);
for i = 1:iterCg
if H1ErrUCg(i) == 0 && i == 1
H1ErrUCg(i) = max(H1ErrUCg);
elseif H1ErrUCg(i) == 0
H1ErrUCg(i) = H1ErrUCg(i-1);
end
end
for i = 1:iterCg
if L2ErrPCg(i) == 0 && i == 1
L2ErrPCg(i) = max(L2ErrPCg);
elseif L2ErrPCg(i) == 0
L2ErrPCg(i) = L2ErrPCg(i-1);
end
end
for i = 1:iterCg
if L2ErrGradUandPCg(i) == 0 && i == 1
L2ErrGradUandPCg(i) = max(L2ErrGradUandPCg);
elseif L2ErrGradUandPCg(i) == 0
L2ErrGradUandPCg(i) = L2ErrGradUandPCg(i-1);
end
end
p = res_conjgrad(size(eqn.A,1)+1:end);
u = res_conjgrad(1:size(eqn.A,1));
time_BPCG_finish = toc;
function [inner_pro] = inner_x_y(C, x,y)
%C = [A 0 ; 0 1]
%(x,y) = (Cx,y) = (Ax_1,y_1) + (x_2,y_2) = transpose(x_1)*A*y_1 +transpose(x_2)*y_2
%trans(y_1)*x_2
inner_pro = transpose(y) * C * x;
end
function [x, i, relres, H1ErrU, L2ErrP, L2ErrGradUandP] = conjgradNew3(~, C,F,~, x, maxiter, tol, u, p, h, M)
z_old = x;
p_old = F - M*z_old;
r_old=p_old;
for i = 1:maxiter
alpha_i = inner_x_y(C, r_old,p_old) / inner_x_y(C, M*p_old,p_old);
z_new = z_old + alpha_i*p_old;
r_new = F - M*z_new;
norm_r = transpose(r_new)*r_new;%inner_x_y(C, r_new,r_new);
H1ErrU(i) = norm(gradient(u - z_new(1:length(u))));
L2ErrP(i) = norm(p - z_new(length(u)+1:end));
L2ErrGradUandP(i) = sqrt((H1ErrU(i))^2 + (L2ErrP(i))^2);
%L2ErrGradUandP(i) = sqrt((norm(gradient(z_new(1:length(uI)),h) - gradient(uI,h)))^2 + (norm(pI - z_new(length(uI)+1:end)))^2);
relres(i) = norm_r;
fprintf('#dof: %8.0u, Bramble Pasciak CG iter: %2.0u, err = %12.8g\n \n',...
length(z_new)+length(z_new(length(u)+1:end)), i, norm_r);
if (norm_r)< tol
break
end
beta_i = inner_x_y(C, M*r_new,p_old) / inner_x_y(C, M*p_old,p_old);
p_new = r_new - beta_i*p_old;
% Update
z_old = z_new;
r_old = r_new;
p_old = p_new;
end
end
function pde = StokesModelMustafaNew
%% Custom Stokes equation model
%
pde = struct('f', 0, 'exactp', @exactp, 'exactu', @exactu,'g_D',@g_D,'g_N',@g_N);
% exact velocity
function z = exactu(p)
x = p(:,1);
y = p(:,2);
z(:,1) = 1/(5*pi()^2).*sin(pi()*x).*sin(2*pi()*y);
z(:,2) = 1/(5*pi()^2).*sin(2*pi()*x).*sin(pi()*y);
end
% exact pressure
function z = exactp(p)
x = p(:,1); y = p(:,2);
z = 2/3-x.^2-y.^2;
end
% Dirichlet boundary condition of velocity
function z = g_D(p)
z = exactu(p);
end
function z = g_N(p)
z(:,1) = 2*ones(size(p,1),1);
z(:,2) = zeros(size(p,1),1);
end
end
Mostafa Kassoum
2024년 1월 12일
The following code is ExtractBoundaryEdges.
function [Dirichlet,bdFlag] = ExtractBoundaryEdges(elem, node)
bdFlag = setboundary(node,elem,'Dirichlet','(x==-1) | (x==1) | (y==1) | (y==-1)');%,'Neumann', '(y==1)');
% bdFlag = setboundary(node,elem,'Neumann','(y==1)');
% [node,elem,bdFlag] = uniformbisect(node,elem,bdFlag);
% [node,elem,bdFlag] = uniformbisect(node,elem,bdFlag);
allEdge = [elem(:,[2,3]); elem(:,[3,1]); elem(:,[1,2])];
Dirichlet = allEdge((bdFlag(:) == 1),:);
end
Walter Roberson
2024년 1월 12일
Need StokesisoP2P0
%bring in external tools
sqm_url = "https://raw.githubusercontent.com/lyc102/ifem/master/mesh/squaremesh.m";
shm_url = "https://raw.githubusercontent.com/lyc102/ifem/master/tool/showmesh.m";
unir_url = "https://raw.githubusercontent.com/lyc102/ifem/master/mesh/uniformrefine.m";
myu_url = "https://raw.githubusercontent.com/lyc102/ifem/master/tool/myunique.m";
setb_url = "https://raw.githubusercontent.com/lyc102/ifem/master/fem/setboundary.m";
tn = tempname();
sqm_localname = fullfile(tn, "squaremesh.m");
shm_localname = fullfile(tn, "showmesh.m");
unir_localname = fullfile(tn, "uniformrefine.m");
myu_localname = fullfile(tn, "myunique.m");
setb_localname = fullfile(tn, "setboundary.m");
mkdir(tn);
addpath(tn)
websave(sqm_localname, sqm_url);
websave(shm_localname, shm_url);
websave(unir_localname, unir_url);
websave(myu_localname, myu_url);
websave(setb_localname, setb_url);
%done bringing in external tools
%Define problem parameters
mu = 1.0; % viscosity
%f = [1; 0; 1; 2]; % forcing term for velocity
%g = [2; 1]; % forcing term for pressure
%% Mesh creation
h = 0.5; % discretization step
[node,elem] = squaremesh([-1,1,-1,1],h);
subplot(1,2,1); showmesh(node,elem);
N = size(node,1);
NTc = size(elem,1);
% Uniform refinement to get a fine grid
[node,elem] = uniformrefine(node,elem);
% Uniform refinement to get a fine grid
subplot(1,2,2); showmesh(node,elem);
%%extract boundary conditions
[Dirichlet, bdFlag] = ExtractBoundaryEdges(elem, node);
%% Load Stokes Data Partial Differential Equations
%Different velocity and pressure equation models
pde = StokesModelMustafaNew;%StokesModelMustafa;%StokesModelMustafaNew;%Stokesdata3;%StokesModelMustafa2;%StokesModelMustafaSin;%;
option.elemType = 'isoP2P0';
option.fquadorder = 3;
tolerance = 1e-3;
max_iterations = 1000;
mesh.node = node;
mesh.elem = elem;
mesh.bdFlag = bdFlag;
%The following line of code is only run to get the eqn (equation structure) easier without defining it every time
[soln,eqn,info] = StokesisoP2P0(node,elem,bdFlag,pde, option);
Unrecognized function or variable 'StokesisoP2P0'.
%% Call the solvers and calulate L2 errors, run time and outcome as well as convergence rate
% [soln,eqn,err,time,solver] = femStokes(mesh,pde,option);
% return
uI = pde.exactu([node; (node(eqn.edge(:,1),:)+node(eqn.edge(:,2),:))/2]);
uI = reshape(uI, [2*length(uI),1]);
pI = pde.exactp([node(elem(:,1),1) node(elem(:,3),2)]);
% Define matrix H
%H = [1 0 0 2; 0 2 1 0; 2 3 1 0; 0 1 0 1];
% Define matrix A
%A = H * H';
% Define matrix B
%B = [1 2 1 0; 0 1 0 0];
% Assemble the finite element system
%[A, B, ~, ~] = assempde(pde,node,elem,'f',f,'g',g,'mu',mu);
tic;
% Assemble the Stokes matrix
N = size(eqn.A,1);
M = size(eqn.B,1);
%StokesMatrix = [eqn.A, eqn.B'; eqn.B, sparse(M,M)];
% Define dimensions of eqn.B'
eqn.B_transpose = eqn.B'; % Transpose of eqn.B
% Construct a sparse MxM matrix (adjust M to your desired size)
M = size(eqn.B, 1); % Number of columns in eqn.B (or rows in eqn.B')
sparse_M_M = sparse(M, M);
% Check dimensions of matrices involved in the problematic line
size_eqn.B = size(eqn.B);
size_eqn.A = size(eqn.A);
% Perform the matrix multiplication explicitly and verify dimensions
result = 2 * eqn.B / (eqn.A)* eqn.B_transpose;
size_result = size(result);
% Concatenate eqn.A, eqn.B_transpose, eqn.B, and sparse_M_M
StokesMatrix = [eqn.A, eqn.B_transpose; eqn.B, sparse_M_M];
% Define matrix M
%M = [2*eye(size(eqn.A,1)) 2*(eqn.A)\eqn.B_transpose; eqn.B 2*eqn.B*(eqn.A)\eqn.B_transpose];
M_top = [2 * eye(size(eqn.A, 1)), 2 * (eqn.A\eqn.B_transpose)];
M_bottom = [eqn.B, 2 * eqn.B * (eqn.A \ eqn.B_transpose)];
M = [M_top; M_bottom];
% Define matrix F
F = [2 * eqn.A \ eqn.f; 2 * eqn.B * (eqn.A \ eqn.f) - eqn.g];
% Define the right-hand side
rhs = [eqn.f; eqn.g];
x0 = ones(size(F));
% Solve the Stokes system
solution = StokesMatrix \ rhs;
sol = M \ F;
% Extract the velocity and pressure components from the solution
u = solution(1:N);
p = solution(N+1:end);
% Calculate the residual R
R = F - M * sol;
% Define matrix C
C = [0.5 * eqn.A zeros(size(transpose(eqn.B))); zeros(size(eqn.B)) eye(size(eqn.B,1))];
% Define an appropriate non-zero vector x
[res_conjgrad, iterCg, relresCg, H1ErrUCg, L2ErrPCg, L2ErrGradUandPCg] = conjgradNew3(eqn.A, C,F ,[eqn.f; eqn.g],x0 , max_iterations, tolerance, u, p, h, M);
for i = 1:iterCg
if H1ErrUCg(i) == 0 && i == 1
H1ErrUCg(i) = max(H1ErrUCg);
elseif H1ErrUCg(i) == 0
H1ErrUCg(i) = H1ErrUCg(i-1);
end
end
for i = 1:iterCg
if L2ErrPCg(i) == 0 && i == 1
L2ErrPCg(i) = max(L2ErrPCg);
elseif L2ErrPCg(i) == 0
L2ErrPCg(i) = L2ErrPCg(i-1);
end
end
for i = 1:iterCg
if L2ErrGradUandPCg(i) == 0 && i == 1
L2ErrGradUandPCg(i) = max(L2ErrGradUandPCg);
elseif L2ErrGradUandPCg(i) == 0
L2ErrGradUandPCg(i) = L2ErrGradUandPCg(i-1);
end
end
p = res_conjgrad(size(eqn.A,1)+1:end);
u = res_conjgrad(1:size(eqn.A,1));
time_BPCG_finish = toc;
function [inner_pro] = inner_x_y(C, x,y)
%C = [A 0 ; 0 1]
%(x,y) = (Cx,y) = (Ax_1,y_1) + (x_2,y_2) = transpose(x_1)*A*y_1 +transpose(x_2)*y_2
%trans(y_1)*x_2
inner_pro = transpose(y) * C * x;
end
function [x, i, relres, H1ErrU, L2ErrP, L2ErrGradUandP] = conjgradNew3(~, C,F,~, x, maxiter, tol, u, p, h, M)
z_old = x;
p_old = F - M*z_old;
r_old=p_old;
for i = 1:maxiter
alpha_i = inner_x_y(C, r_old,p_old) / inner_x_y(C, M*p_old,p_old);
z_new = z_old + alpha_i*p_old;
r_new = F - M*z_new;
norm_r = transpose(r_new)*r_new;%inner_x_y(C, r_new,r_new);
H1ErrU(i) = norm(gradient(u - z_new(1:length(u))));
L2ErrP(i) = norm(p - z_new(length(u)+1:end));
L2ErrGradUandP(i) = sqrt((H1ErrU(i))^2 + (L2ErrP(i))^2);
%L2ErrGradUandP(i) = sqrt((norm(gradient(z_new(1:length(uI)),h) - gradient(uI,h)))^2 + (norm(pI - z_new(length(uI)+1:end)))^2);
relres(i) = norm_r;
fprintf('#dof: %8.0u, Bramble Pasciak CG iter: %2.0u, err = %12.8g\n \n',...
length(z_new)+length(z_new(length(u)+1:end)), i, norm_r);
if (norm_r)< tol
break
end
beta_i = inner_x_y(C, M*r_new,p_old) / inner_x_y(C, M*p_old,p_old);
p_new = r_new - beta_i*p_old;
% Update
z_old = z_new;
r_old = r_new;
p_old = p_new;
end
end
function pde = StokesModelMustafaNew
%% Custom Stokes equation model
%
pde = struct('f', 0, 'exactp', @exactp, 'exactu', @exactu,'g_D',@g_D,'g_N',@g_N);
% exact velocity
function z = exactu(p)
x = p(:,1);
y = p(:,2);
z(:,1) = 1/(5*pi()^2).*sin(pi()*x).*sin(2*pi()*y);
z(:,2) = 1/(5*pi()^2).*sin(2*pi()*x).*sin(pi()*y);
end
% exact pressure
function z = exactp(p)
x = p(:,1); y = p(:,2);
z = 2/3-x.^2-y.^2;
end
% Dirichlet boundary condition of velocity
function z = g_D(p)
z = exactu(p);
end
function z = g_N(p)
z(:,1) = 2*ones(size(p,1),1);
z(:,2) = zeros(size(p,1),1);
end
end
function [Dirichlet,bdFlag] = ExtractBoundaryEdges(elem, node)
bdFlag = setboundary(node,elem,'Dirichlet','(x==-1) | (x==1) | (y==1) | (y==-1)');%,'Neumann', '(y==1)');
% bdFlag = setboundary(node,elem,'Neumann','(y==1)');
% [node,elem,bdFlag] = uniformbisect(node,elem,bdFlag);
% [node,elem,bdFlag] = uniformbisect(node,elem,bdFlag);
allEdge = [elem(:,[2,3]); elem(:,[3,1]); elem(:,[1,2])];
Dirichlet = allEdge((bdFlag(:) == 1),:);
end
Mostafa Kassoum
2024년 1월 13일
You need download the codes in the following link on MATLAB to get what you need.
Walter Roberson
2024년 1월 13일
%bring in external tools
sqm_url = "https://raw.githubusercontent.com/lyc102/ifem/master/mesh/squaremesh.m";
shm_url = "https://raw.githubusercontent.com/lyc102/ifem/master/tool/showmesh.m";
unir_url = "https://raw.githubusercontent.com/lyc102/ifem/master/mesh/uniformrefine.m";
myu_url = "https://raw.githubusercontent.com/lyc102/ifem/master/tool/myunique.m";
setb_url = "https://raw.githubusercontent.com/lyc102/ifem/master/fem/setboundary.m";
S2P0_url = "https://raw.githubusercontent.com/lyc102/ifem/master/equation/StokesisoP2P0.m";
dof_url = "https://raw.githubusercontent.com/lyc102/ifem/master/dof/dof3P2.m";
tn = tempname();
sqm_localname = fullfile(tn, "squaremesh.m");
shm_localname = fullfile(tn, "showmesh.m");
unir_localname = fullfile(tn, "uniformrefine.m");
myu_localname = fullfile(tn, "myunique.m");
setb_localname = fullfile(tn, "setboundary.m");
S2P0_localname = fullfile(tn, "StokesisoP2P0.m");
dof3_localname = fullfile(tn, "dof3P2.m");
dof_localname = fullfile(tn, "dofP2.m");
mkdir(tn);
addpath(tn)
websave(sqm_localname, sqm_url);
websave(shm_localname, shm_url);
websave(unir_localname, unir_url);
websave(myu_localname, myu_url);
websave(setb_localname, setb_url);
websave(S2P0_localname, S2P0_url);
websave(dof_localname, dof_url);
websave(dof3_localname, dof_url);
%done bringing in external tools
%Define problem parameters
mu = 1.0; % viscosity
%f = [1; 0; 1; 2]; % forcing term for velocity
%g = [2; 1]; % forcing term for pressure
%% Mesh creation
h = 0.5; % discretization step
[node,elem] = squaremesh([-1,1,-1,1],h);
subplot(1,2,1); showmesh(node,elem);
N = size(node,1);
NTc = size(elem,1);
% Uniform refinement to get a fine grid
[node,elem] = uniformrefine(node,elem);
% Uniform refinement to get a fine grid
subplot(1,2,2); showmesh(node,elem);
%%extract boundary conditions
[Dirichlet, bdFlag] = ExtractBoundaryEdges(elem, node);
%% Load Stokes Data Partial Differential Equations
%Different velocity and pressure equation models
pde = StokesModelMustafaNew;%StokesModelMustafa;%StokesModelMustafaNew;%Stokesdata3;%StokesModelMustafa2;%StokesModelMustafaSin;%;
option.elemType = 'isoP2P0';
option.fquadorder = 3;
tolerance = 1e-3;
max_iterations = 1000;
mesh.node = node;
mesh.elem = elem;
mesh.bdFlag = bdFlag;
%The following line of code is only run to get the eqn (equation structure) easier without defining it every time
[soln,eqn,info] = StokesisoP2P0(node,elem,bdFlag,pde, option);
Index in position 2 exceeds array bounds. Index must not exceed 3.
Error in dofP2 (line 17)
totalEdge = uint32([elem(:,[1 2]); elem(:,[1 3]); elem(:,[1 4]); ...
Error in StokesisoP2P0 (line 25)
[tempvar,edgeC] = dofP2(elemC);
%% Call the solvers and calulate L2 errors, run time and outcome as well as convergence rate
% [soln,eqn,err,time,solver] = femStokes(mesh,pde,option);
% return
uI = pde.exactu([node; (node(eqn.edge(:,1),:)+node(eqn.edge(:,2),:))/2]);
uI = reshape(uI, [2*length(uI),1]);
pI = pde.exactp([node(elem(:,1),1) node(elem(:,3),2)]);
% Define matrix H
%H = [1 0 0 2; 0 2 1 0; 2 3 1 0; 0 1 0 1];
% Define matrix A
%A = H * H';
% Define matrix B
%B = [1 2 1 0; 0 1 0 0];
% Assemble the finite element system
%[A, B, ~, ~] = assempde(pde,node,elem,'f',f,'g',g,'mu',mu);
tic;
% Assemble the Stokes matrix
N = size(eqn.A,1);
M = size(eqn.B,1);
%StokesMatrix = [eqn.A, eqn.B'; eqn.B, sparse(M,M)];
% Define dimensions of eqn.B'
eqn.B_transpose = eqn.B'; % Transpose of eqn.B
% Construct a sparse MxM matrix (adjust M to your desired size)
M = size(eqn.B, 1); % Number of columns in eqn.B (or rows in eqn.B')
sparse_M_M = sparse(M, M);
% Check dimensions of matrices involved in the problematic line
size_eqn.B = size(eqn.B);
size_eqn.A = size(eqn.A);
% Perform the matrix multiplication explicitly and verify dimensions
result = 2 * eqn.B / (eqn.A)* eqn.B_transpose;
size_result = size(result);
% Concatenate eqn.A, eqn.B_transpose, eqn.B, and sparse_M_M
StokesMatrix = [eqn.A, eqn.B_transpose; eqn.B, sparse_M_M];
% Define matrix M
%M = [2*eye(size(eqn.A,1)) 2*(eqn.A)\eqn.B_transpose; eqn.B 2*eqn.B*(eqn.A)\eqn.B_transpose];
M_top = [2 * eye(size(eqn.A, 1)), 2 * (eqn.A\eqn.B_transpose)];
M_bottom = [eqn.B, 2 * eqn.B * (eqn.A \ eqn.B_transpose)];
M = [M_top; M_bottom];
% Define matrix F
F = [2 * eqn.A \ eqn.f; 2 * eqn.B * (eqn.A \ eqn.f) - eqn.g];
% Define the right-hand side
rhs = [eqn.f; eqn.g];
x0 = ones(size(F));
% Solve the Stokes system
solution = StokesMatrix \ rhs;
sol = M \ F;
% Extract the velocity and pressure components from the solution
u = solution(1:N);
p = solution(N+1:end);
% Calculate the residual R
R = F - M * sol;
% Define matrix C
C = [0.5 * eqn.A zeros(size(transpose(eqn.B))); zeros(size(eqn.B)) eye(size(eqn.B,1))];
% Define an appropriate non-zero vector x
[res_conjgrad, iterCg, relresCg, H1ErrUCg, L2ErrPCg, L2ErrGradUandPCg] = conjgradNew3(eqn.A, C,F ,[eqn.f; eqn.g],x0 , max_iterations, tolerance, u, p, h, M);
for i = 1:iterCg
if H1ErrUCg(i) == 0 && i == 1
H1ErrUCg(i) = max(H1ErrUCg);
elseif H1ErrUCg(i) == 0
H1ErrUCg(i) = H1ErrUCg(i-1);
end
end
for i = 1:iterCg
if L2ErrPCg(i) == 0 && i == 1
L2ErrPCg(i) = max(L2ErrPCg);
elseif L2ErrPCg(i) == 0
L2ErrPCg(i) = L2ErrPCg(i-1);
end
end
for i = 1:iterCg
if L2ErrGradUandPCg(i) == 0 && i == 1
L2ErrGradUandPCg(i) = max(L2ErrGradUandPCg);
elseif L2ErrGradUandPCg(i) == 0
L2ErrGradUandPCg(i) = L2ErrGradUandPCg(i-1);
end
end
p = res_conjgrad(size(eqn.A,1)+1:end);
u = res_conjgrad(1:size(eqn.A,1));
time_BPCG_finish = toc;
function [inner_pro] = inner_x_y(C, x,y)
%C = [A 0 ; 0 1]
%(x,y) = (Cx,y) = (Ax_1,y_1) + (x_2,y_2) = transpose(x_1)*A*y_1 +transpose(x_2)*y_2
%trans(y_1)*x_2
inner_pro = transpose(y) * C * x;
end
function [x, i, relres, H1ErrU, L2ErrP, L2ErrGradUandP] = conjgradNew3(~, C,F,~, x, maxiter, tol, u, p, h, M)
z_old = x;
p_old = F - M*z_old;
r_old=p_old;
for i = 1:maxiter
alpha_i = inner_x_y(C, r_old,p_old) / inner_x_y(C, M*p_old,p_old);
z_new = z_old + alpha_i*p_old;
r_new = F - M*z_new;
norm_r = transpose(r_new)*r_new;%inner_x_y(C, r_new,r_new);
H1ErrU(i) = norm(gradient(u - z_new(1:length(u))));
L2ErrP(i) = norm(p - z_new(length(u)+1:end));
L2ErrGradUandP(i) = sqrt((H1ErrU(i))^2 + (L2ErrP(i))^2);
%L2ErrGradUandP(i) = sqrt((norm(gradient(z_new(1:length(uI)),h) - gradient(uI,h)))^2 + (norm(pI - z_new(length(uI)+1:end)))^2);
relres(i) = norm_r;
fprintf('#dof: %8.0u, Bramble Pasciak CG iter: %2.0u, err = %12.8g\n \n',...
length(z_new)+length(z_new(length(u)+1:end)), i, norm_r);
if (norm_r)< tol
break
end
beta_i = inner_x_y(C, M*r_new,p_old) / inner_x_y(C, M*p_old,p_old);
p_new = r_new - beta_i*p_old;
% Update
z_old = z_new;
r_old = r_new;
p_old = p_new;
end
end
function pde = StokesModelMustafaNew
%% Custom Stokes equation model
%
pde = struct('f', 0, 'exactp', @exactp, 'exactu', @exactu,'g_D',@g_D,'g_N',@g_N);
% exact velocity
function z = exactu(p)
x = p(:,1);
y = p(:,2);
z(:,1) = 1/(5*pi()^2).*sin(pi()*x).*sin(2*pi()*y);
z(:,2) = 1/(5*pi()^2).*sin(2*pi()*x).*sin(pi()*y);
end
% exact pressure
function z = exactp(p)
x = p(:,1); y = p(:,2);
z = 2/3-x.^2-y.^2;
end
% Dirichlet boundary condition of velocity
function z = g_D(p)
z = exactu(p);
end
function z = g_N(p)
z(:,1) = 2*ones(size(p,1),1);
z(:,2) = zeros(size(p,1),1);
end
end
function [Dirichlet,bdFlag] = ExtractBoundaryEdges(elem, node)
bdFlag = setboundary(node,elem,'Dirichlet','(x==-1) | (x==1) | (y==1) | (y==-1)');%,'Neumann', '(y==1)');
% bdFlag = setboundary(node,elem,'Neumann','(y==1)');
% [node,elem,bdFlag] = uniformbisect(node,elem,bdFlag);
% [node,elem,bdFlag] = uniformbisect(node,elem,bdFlag);
allEdge = [elem(:,[2,3]); elem(:,[3,1]); elem(:,[1,2])];
Dirichlet = allEdge((bdFlag(:) == 1),:);
end
Mostafa Kassoum
2024년 1월 13일
It worked for me anf I got the results for everything and they are nice results except the error of pressure P, they are so big. I think there is a mistake but I didn't got it.
Mostafa Kassoum
2024년 1월 13일
Foe example I got the following results:
H1 error for velocity u are: 0.129922159092270 0.138027304457401
L2 error for pressure p are: 12.6454430756918 12.6450517755845
I need know why are L2 error p big
Walter Roberson
2024년 1월 13일
You can copy my entire post to your system and test it locally.
I suspect that something you have installed does not match the latest code at https://github.com/lyc102/ifem
Mostafa Kassoum
2024년 1월 23일
I copied your entire post and runned it but I got also a big value of L2 error 11.3255614294248 11.3251245255107
As you can see they are almost equal to the values I obtained previously. Where do you think the mistake is?
답변 (0개)
참고 항목
카테고리
Help Center 및 File Exchange에서 Geometry and Mesh에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!오류 발생
페이지가 변경되었기 때문에 동작을 완료할 수 없습니다. 업데이트된 상태를 보려면 페이지를 다시 불러오십시오.
웹사이트 선택
번역된 콘텐츠를 보고 지역별 이벤트와 혜택을 살펴보려면 웹사이트를 선택하십시오. 현재 계신 지역에 따라 다음 웹사이트를 권장합니다:
또한 다음 목록에서 웹사이트를 선택하실 수도 있습니다.
사이트 성능 최적화 방법
최고의 사이트 성능을 위해 중국 사이트(중국어 또는 영어)를 선택하십시오. 현재 계신 지역에서는 다른 국가의 MathWorks 사이트 방문이 최적화되지 않았습니다.
미주
- América Latina (Español)
- Canada (English)
- United States (English)
유럽
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)
아시아 태평양
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)