アダマール変換の定義 F(u,v) &= \frac{1}{n} \sum_{x=0}^{n-1}\sum_{y=0}^{n-1} f(x,y)\cdot (-1)^{q(x,y,u,v)} \\ q(x,y,u,v) &= \sum_{i=0}^{b-1}x_i g_i(u) + y_i g_i (v)\\ g_i(u) &= u_{b-i} + u_{b-i-1} ただしbはu(およびv)を２進数で表現したときに…

DCT-IIの定義式 ただし function show_dct_basis clear; close all; m = 8; N = m*m; function F=getBasisF(u,v) F = zeros(m,m); for y=0:m-1 for x=0:m-1 F(y+1,x+1) = sqrt(2/m) * calcC(u)*calcC(v) ... * cos(((2*x+1)*u*pi)/(2*m)) ... * cos(((2*y+1)…

function ycbcr_test close all; clear; %% lena im = imread('lena.jpg'); [X,Y,c] = size(im); Ys= zeros(X,Y); Cbs= zeros(X,Y); Crs = zeros(X,Y); for y=1:Y for x=1:X R = im(y,x,1); G = im(y,x,1); B = im(y,x,1); [Y,Cb,Cr] = RGB2YCbCR(R,G,B); Ys…

function encoder close all; clear; %% initialization N = 256; im = zeros(N,N,3); im2 = zeros(N,N,3); Ys = zeros(N,N); Cbs = zeros(N,N); Crs = zeros(N,N); %% Main routine im = makeRGB(im); subplot(141); imshow(im); title("original RGB"); %i…

緑が強い部分が白くなる function ycbcr_test close all; clear; %% lena im = imread('lena.jpg'); [X,Y,c] = size(im); im2 = zeros(X,Y,c); for y=1:Y for x=1:X R = im(y,x,1); G = im(y,x,1); B = im(y,x,1); [Y,Cb,Cr] = RGB2YCbCR(R,G,B); im2(y,x,1)…

function chebyshev figure();hold on; title('Chebyshev polynomial'); for n=1:4 theta = linspace(0,2*pi,100); x = cos(theta); Tnx = cos(n*theta); plot(x,Tnx); end l = legend({'N=1' 'N=2' 'N=3' 'N=4'}); set(l,'FontSize',14); end チェビシェフ…

結構めんどくさい うまく動いてない途中で三重対角行列を係数行列とする連立１次方程式の解を求める必要があったが、TDMA（三重対角行列アルゴリズム）という解法がある。 function spline_interpolation clear all; close all; format compact; div = 100; …

N=3ですでにサインカーブがうまく補間されている function lagrange_interpolation clear all; close all; div = 100; x = linspace(0,pi,100); y = sin(x); plot(x,y); hold on; a = x(1); b = x(end); fx = (x-b)/(a-b)*sin(a) + (x-a)/(b-a)*sin(b); plot…

http://www.takuichi.net/hobby/edu/em/mom/numerical/index-j.htm を計算する リーマン積分（区分求積法） function riemann clear all; close all; divlist = [10e2,10e3,10e4,10e5,10e6,10e7]; errorlist = []; for cycle=1:6 div = divlist(cycle); x = …

function conway2d clear all; close all; CAsize = 100; [X,Y] = meshgrid(1,CAsize+1); figure(); hold on; plot(X,Y,'k'); plot(Y,X,'k'); I = zeros(CAsize+1); surface(I); colormap(gray); axis off; I = rand(CAsize+1)>0.85; I(1,:) = 0; I(CAsize,:…

1-256の規則を全て適用してみた function wolfram1d_2 clear all; close all; CAsize = 100; [X,Y] = meshgrid(1,CAsize+1); k = {'111','110','101','100','011','010','001','000'}; % v = {'1','0','1','0','0','1','0','1'}; %165 % writerObj = VideoWr…

規則182を適用 function wolfram1d_2 clear all; close all; CAsize = 100; [X,Y] = meshgrid(1,CAsize+1); figure(); hold on; plot(X,Y,'k'); plot(Y,X,'k'); I = zeros(CAsize+1); surface(I); colormap(gray); axis off; origin = zeros(1,CAsize); % fo…

https://www.mathworks.com/help/matlab/visualize/representing-a-matrix-as-a-surface.html figure k = 5; n = 2^k-1; theta = pi*(-n:2:n)/n; phi = (pi/2)*(-n:2:n)'/n; X = cos(phi)*cos(theta); Y = cos(phi)*sin(theta); Z = sin(phi)*ones(size(thet…

function velocity_variable_diffusion clear all; close all; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% CAPTURE = false; VIEWTOP = true; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %condition …

function diff2d_constD clear all; close all; % use surf() function for visualization xmax=2; ymax=2; deltaT = 0.1; % D=0.05; % dx=D*(deltaT)/() % dy= xwidth=31; ywidth=31; x=linspace(-xmax,xmax,xwidth); diffx=diff(x); dx=diffx(1); y=linspa…

function diff2d_constD clear all; close all; % use surf() function for visualization %% define capture CAPTURE = false; xmax=2; ymax=2; deltaT = 0.1; D=0.05; % dx=D*(deltaT)/() % dy= xwidth=31; ywidth=31; x=linspace(-xmax,xmax,xwidth); dif…

function diff1d clear all; close all; width=1; A=1; xmax=5; partition=101; x=linspace(-xmax,xmax,partition); D=0.04; dxs=diff(x); dx=dxs(1); deltaT=0.1; d=D*deltaT/dx^2; tmax=4; global f1; global f2; global f3; global f4; set(gcf, 'Positio…

function advection_x clear all; close all; format shortG; %%define capture CAPTURE = false; partition=101; xmax=20; xx=linspace(0,xmax,partition)'; dxs=diff(xx); dx=dxs(1); deltaT=0.02; u=2; c=u*deltaT/dx; global d; global h; d=0.4; h=1.2;…

function potential_doublet clear all; close all; x=-2:0.1:2; y=-2:0.1:2; [X,Y] = meshgrid(x,y); k=1.0; Phi = k*X./(X.^2+Y.^2); Psi = -k*Y./(X.^2+Y.^2); contour(X,Y,Phi,50,'LineColor','red'); hold on; contour(X,Y,Psi,50,'LineColor','blue');…

function uniform_doublet_vortex clear all; close all; x=-2.0:0.1:2.0; y=-2.0:0.1:2.0; [X,Y] = meshgrid(x,y); U=2.0; alpha_deg=0; alpha=alpha_deg/180.0*pi; R=0.8; r=sqrt(X.^2+Y.^2); theta = atan2(Y,X); gamma = - 5.0; Phi = U*(r+R^2./r).*cos…

計算条件 速度ポテンシャルを求める function laplace2d_velocity_potential clear all; close all; format shortG; %%define capture CAPTURE = true; global L;global U;global W;global Xobs;global WobsX;global WobsY; L = 2.0; U = 1.0; W = 1.8; Xobs…

計算条件 流れ関数の計算 function laplace2d clear all; close all; %% define capture CAPTURE = false; Xobs = 0.8; WobsX = 0.4; WobsY = 0.5; L = 2; W = 1; deltaX = 100; deltaY = deltaX*W/2; x=linspace(0,L,deltaX); y=linspace(0,W,deltaY); [X,Y…

function wing_flow_conformal clear all;close all; R=0.8; U=10; alpha_deg = 20; alpha = alpha_deg/180*pi; gamma = -4*pi*R*U*sin(alpha); rho = 10; x=-2:0.05:2; y=-2:0.05:2; [X,Y] = meshgrid(x,y); Xi= X +i*Y; xi0 = -0.12; eta0 = 0.52; Rdash =…

function plate_flow_conformal clear all;close all; R=0.8; U=10; alpha_deg = 20; alpha = alpha_deg/180*pi; gamma = -4*pi*R*U*sin(alpha); rho = 10; x=-3:0.05:3; y=-3:0.05:3; [X,Y] = meshgrid(x,y); Xi= X +i*Y; % W = U*(Xi + R^2./Xi*exp(-2*i*a…

function pressure_map clear all; close all; def = 2; x=-def:0.1:def; y=-def:0.1:def; [X,Y] = meshgrid(x,y); U= 2.0; alpha_deg=0; alpha=alpha_deg/180.0*pi; R=0.8; r=sqrt(X.^2+Y.^2); theta = atan2(Y,X); gamma = -15.0; rho = 10.0; Phi = U*(r+…

function doublet_in_uniform_flow % clear all; % close all; x=-2.0:0.1:2.0; y=-2.0:0.1:2.0; [X,Y] = meshgrid(x,y); U=1.0; alpha_deg=0; alpha=alpha_deg/180.0*pi; R=0.8; r=sqrt(X.^2+Y.^2); theta = atan2(Y,X); gamma = 1.0; Phi = U*(r+R^2./r).*…

function doublet_in_uniform_flow % clear all; % close all; x=-2.0:0.1:2.0; y=-2.0:0.1:2.0; [X,Y] = meshgrid(x,y); U=1.0; alpha_deg=0; alpha=alpha_deg/180.0*pi; R=0.8; r=sqrt(X.^2+Y.^2); theta = atan2(Y,X); gamma = 1.0; Phi = U*(r+R^2./r).*…

function potential_source close all; clear all; x=-10:0.1:10; y=-10:0.1:10; [X,Y] = meshgrid(x,y); Q=1.0; r=sqrt(X.^2+Y.^2); U=1.0; alpha_deg=5; alpha=alpha_deg/180.0*pi; Phi = Q/(2*pi) * log(r); Psi = Q/(2*pi) * atan2(Y,X); contour(Phi,20…

function potential_source close all; clear all; x=-10:0.1:10; y=-10:0.1:10; [X,Y] = meshgrid(x,y); Q=1.0; r=sqrt(X.^2+Y.^2); U=1.0; alpha_deg=5; alpha=alpha_deg/180.0*pi; Phi = Q/(2*pi) * log(r); Psi = Q/(2*pi) * atan2(Y,X); contour(Phi,20…

function potential_flow x=0:0.1:10; y=0:0.1:10; [X,Y] = meshgrid(x,y); U=1.0; alpha_deg=5; alpha=alpha_deg/180.0*pi; Phi = U*(X.*cos(alpha)+Y.*sin(alpha)); Psi = U*(Y.*cos(alpha)-X.*sin(alpha)); contour(Phi,'LineColor','red'); hold on; con…