日韩av黄I国产麻豆传媒I国产91av视频在线观看I日韩一区二区三区在线看I美女国产在线I麻豆视频国产在线观看I成人黄色短片

歡迎訪問 生活随笔!

生活随笔

當(dāng)前位置: 首頁(yè) >

利用matlab绘制二维均匀流线和向量场

發(fā)布時(shí)間:2023/12/20 46 豆豆
生活随笔 收集整理的這篇文章主要介紹了 利用matlab绘制二维均匀流线和向量场 小編覺得挺不錯(cuò)的,現(xiàn)在分享給大家,幫大家做個(gè)參考.

利用matlab繪制二維均勻流線和向量場(chǎng)(向量場(chǎng)彩色箭頭,顏色隨變量變化)

  • 0前言
  • 1 均勻流線的繪制
  • 2 繪制彩色的短線圖
  • 3 繪制彩色的均勻流線
  • 4 運(yùn)動(dòng)的彩色箭頭流線圖

0前言

之前一篇文章matlab流場(chǎng)可視化后處理我簡(jiǎn)單介紹了很多繪制流場(chǎng)可視化的方法,但是在流線方面一直不是很滿意,因?yàn)閙atlab自帶的繪制流線函數(shù)streamslice沒有開源,有些參數(shù)(比如每根流線的顏色)不能修改等,所以我嘗試著自己寫一種繪制二維均勻流線的程序,以便滿足更多的需求。

其中,matlab官方file exchange里,也有一些大神的代碼給出了很好的結(jié)果:
Evenly Spaced Streamlines

Streamcolor

StreakArrow

關(guān)于如何用matlab繪制均勻的流線,我采用的方法是下面這篇論文提到的算法:
Jobard B., Lefer W. (1997) Creating Evenly-Spaced Streamlines of Arbitrary Density. In: Lefer W., Grave M. (eds) Visualization in Scientific Computing ’97. Eurographics. Springer, Vienna

1 均勻流線的繪制

論文Creating Evenly-Spaced Streamlines of Arbitrary Density中提出的繪制流線的算法可以描述為:

先初始化一根流線,然后在這個(gè)流線周圍逐漸的擴(kuò)展新的流線。文章中定義了兩個(gè)距離,一個(gè)是用于產(chǎn)生新流線種子點(diǎn)的距離dstart,一個(gè)是用于結(jié)束流線生成的距離dend。新的流線擴(kuò)展方式為,在距離舊的流線為dstart的隨機(jī)位置上,生成一個(gè)種子點(diǎn),對(duì)該種子點(diǎn)進(jìn)行新流線的生成(同時(shí)向頭和尾段生成流線),當(dāng)新流線兩端的點(diǎn)距某個(gè)舊流線的距離小于dend時(shí),停止生成流線,將該流線放入舊流線內(nèi)保存。當(dāng)任意位置都不能滿足新流線的生成時(shí),則認(rèn)為所有流線都已經(jīng)生成,結(jié)束計(jì)算。

用偽代碼的形式可以表述為

1 隨機(jī)繪制第一根流線,作為初始流線。把初始流線加入已生成的流線列表內(nèi)。把初始流線記為當(dāng)前流線。 2 令Finish=False,當(dāng)Finish==False時(shí),進(jìn)行循環(huán)3 選擇一個(gè)種子點(diǎn),該點(diǎn)距 當(dāng)前流線 的距離為dstart,且距離其它 已生成的流線列表內(nèi) 的流線距離大于dstart4 如果能夠找到一個(gè)符合條件的種子點(diǎn),則生成新的流線。5 新流線以種子點(diǎn)為開始,向前向后伸展。當(dāng)流線的前端新計(jì)算的流線點(diǎn)距離某條流線的距離小于dend,則停止前端流線的生成。后端同樣。6 將新流線保存到已生成的流線列表內(nèi)。7 如果不能找到符合條件的種子點(diǎn),則替換 當(dāng)前流線 ,重復(fù)48 如果 已生成的流線列表內(nèi) 所有流線都經(jīng)過步驟7嘗試一遍,則說明計(jì)算結(jié)束。令Finish==True, 停止2的循環(huán)。

由于在計(jì)算中,用到了大量的距離計(jì)算,比如生成新種子點(diǎn)時(shí),需要不停的和其它所有流線計(jì)算距離,保證大于dstart,又比如在計(jì)算新流線時(shí),新流線每次生成一個(gè)點(diǎn),該點(diǎn)都需要和所有以有的流線進(jìn)行距離計(jì)算,保證距離大于dend。

這種頻繁的距離計(jì)算數(shù)量,在量級(jí)上為N^2,這導(dǎo)致當(dāng)流線上點(diǎn)比較多或者流線比較多時(shí),計(jì)算會(huì)非常緩慢。我們可以利用一些其它的方法減少這種運(yùn)算,比如利用劃分網(wǎng)格的方式,使得每次計(jì)算距離,只需要計(jì)算周圍網(wǎng)格(比如2維的周圍網(wǎng)格有8個(gè))內(nèi)所有點(diǎn)的距離,而無需計(jì)算全局所有點(diǎn)的距離?;蛘呃脴湫嗡阉魉惴?#xff0c;可以將計(jì)算量降到NlgN的量級(jí)上。

當(dāng)然,這里我利用一種更簡(jiǎn)單的方法進(jìn)行計(jì)算。利用劃分網(wǎng)格的思想,將新種子點(diǎn)的生成位置限定在網(wǎng)格上。而且,假設(shè)如果兩個(gè)點(diǎn)在同一個(gè)網(wǎng)格內(nèi),說明兩個(gè)點(diǎn)距離小于d,如果不在同一網(wǎng)格內(nèi),兩個(gè)點(diǎn)距離大于d。這種假設(shè)雖然對(duì)計(jì)算結(jié)果不精準(zhǔn),但是仍然是可接受的,而且不需要計(jì)算距離。

用偽代碼的形式可以表述為

1 生成距離dstart的網(wǎng)格。生成距離dend的網(wǎng)格。網(wǎng)格內(nèi)初始化為0 2 當(dāng)start網(wǎng)格內(nèi)所有區(qū)域都為1,則停止循環(huán)3 尋找start網(wǎng)格區(qū)域?yàn)?的區(qū)域,隨機(jī)選擇一個(gè)區(qū)域作為生成種子點(diǎn)。4 以種子點(diǎn)為起始點(diǎn),向前向后繪制流線5 查詢新生成的流線點(diǎn)對(duì)應(yīng)的end網(wǎng)格,記為當(dāng)前網(wǎng)格。6 如果當(dāng)前網(wǎng)格和上一個(gè)點(diǎn)的當(dāng)前網(wǎng)格是同一個(gè),則認(rèn)為當(dāng)前網(wǎng)格沒有移動(dòng),繼續(xù)生成下一個(gè)流線點(diǎn)7 如果當(dāng)前網(wǎng)格移動(dòng),且當(dāng)前end網(wǎng)格為1,則停止生成流線。8 如果當(dāng)前網(wǎng)格移動(dòng),且當(dāng)前end網(wǎng)格為0,則將當(dāng)前網(wǎng)格標(biāo)記為1,繼續(xù)生成下一個(gè)流線點(diǎn)9 將流線保存。流線經(jīng)過的所有start網(wǎng)格標(biāo)記為1。流線經(jīng)過的所有end網(wǎng)格標(biāo)記為1。之后重復(fù)3。

matlab代碼為:

clear clc close all%載入數(shù)據(jù) load wind N=5; xx=x(:,:,N); yy=y(:,:,N); uu=u(:,:,N); vv=v(:,:,N); %由于這里R2018b之前的版本會(huì)出現(xiàn)兼容性問題,所以更改 %xmin=min(x,[],'all');xmax=max(x,[],'all'); %ymin=min(y,[],'all');ymax=max(y,[],'all'); xmin=min(min(min(x)));xmax=max(max(max(x))); ymin=min(min(min(y)));ymax=max(max(max(y)));%流線求解 streamline_sum=my_streamline(xx,yy,uu,vv,0.05);%繪制流線 figure(1) hold on xlim([xmin,xmax]); ylim([ymin,ymax]); xy_ratio = get(gca, 'DataAspectRatio'); xy_ratio = xy_ratio(1:2); xy_lim = axis;for j=1:size(streamline_sum,2)if size(streamline_sum{j},1)<=1continue%忽略異常流線endplot(streamline_sum{j}(:,1),streamline_sum{j}(:,2),'k');%繪制箭頭arrow_direction=streamline_sum{j}(end,:)-streamline_sum{j}(end-1,:);plot_arrow(xy_lim,xy_ratio,streamline_sum{j}(end,:),[0,0,0],0.8,arrow_direction) end hold offfunction streamline_sum=my_streamline(x,y,u,v,dstart) %0處理前設(shè)置 %設(shè)置網(wǎng)格密度(01區(qū)間內(nèi)歸一化的長(zhǎng)度) %dstart=0.05;默認(rèn)0.05 dend=0.5*dstart;%xmin=min(x,[],'all');xmax=max(x,[],'all'); %ymin=min(y,[],'all');ymax=max(y,[],'all'); xmin=min(min(min(x)));xmax=max(max(max(x))); ymin=min(min(min(y)));ymax=max(max(max(y)));%歸一化,將流場(chǎng)縮放為0-1區(qū)間的矩形 xn=(x-xmin)/(xmax-xmin); yn=(y-ymin)/(ymax-ymin); un=u/(xmax-xmin); vn=v/(ymax-ymin);num_start=ceil((0.5-dstart/2)/dstart)*2+1; num_end=ceil((0.5-dend/2)/dend)*2+1;%初始化所有網(wǎng)格點(diǎn),0代表可以放置新點(diǎn),1代表已經(jīng)存在原有的點(diǎn) xy_start=zeros(num_start,num_start); xy_end=zeros(num_end,num_end);%1當(dāng)xy_start內(nèi)還有可放置的新點(diǎn)的位置時(shí),進(jìn)行循環(huán) k=0;%循環(huán)次數(shù),也是流線個(gè)數(shù) while ~all(xy_start,'all')k=k+1;%2隨機(jī)一個(gè)start內(nèi)網(wǎng)格點(diǎn)作為種子點(diǎn)[start_id_y,start_id_x]=find(xy_start==0);randnum=randi(size(start_id_y,1));x_pos_i=id2axis(dstart,start_id_x(randnum,1));y_pos_i=id2axis(dstart,start_id_y(randnum,1));%3繪制流線streamline_i_1 = stream2(xn,yn, un, vn,x_pos_i,y_pos_i,0.2);streamline_i_2 = stream2(xn,yn,-un,-vn,x_pos_i,y_pos_i,0.2);%4以xy_end為標(biāo)準(zhǔn),刪除自相交或間隔太近的點(diǎn)。并順便標(biāo)記xy_end[streamline_i_1,xy_end,xy_start]=delete_self(streamline_i_1{1},xy_end,dend,xy_start,dstart);[streamline_i_2,xy_end,xy_start]=delete_self(streamline_i_2{1},xy_end,dend,xy_start,dstart);%5保存streamline_k=[flipud(streamline_i_2);streamline_i_1(2:end,:)];%新的流線streamline_sum{k}=[xmin+streamline_k(:,1)*(xmax-xmin),ymin+streamline_k(:,2)*(ymax-ymin)];%從歸一化還原 end endfunction xpoint=id2axis(distance,id) %取網(wǎng)格的中點(diǎn) N=ceil((0.5-distance/2)/distance)*2+1;%分割的數(shù)量 min_distance=(1-(N-2)*distance)/2;%兩端最小的距離 if id==1xpoint=min_distance/2; elseif id==Nxpoint=1-min_distance/2; elsexpoint=min_distance+(id-1.5)*distance; end endfunction [sl_i,xy_end,xy_start]=delete_self(sl_i,xy_end,dend,xy_start,dstart) %sl_i streamline流線,兩列N行形式 N=size(sl_i,1); pos_id_last=axis2id(sl_i(1,1),sl_i(1,2),dend); xy_end(pos_id_last)=1;%第一個(gè)點(diǎn)標(biāo)記%順便標(biāo)記xy_start pos_id_s=axis2id(sl_i(1,1),sl_i(1,2),dstart); xy_start(pos_id_s)=1; for j=2:Npos_id_now=axis2id(sl_i(j,1),sl_i(j,2),dend);if pos_id_now~=pos_id_last%如果現(xiàn)在的點(diǎn)和原有的點(diǎn)在同一區(qū)域,則不管它%如果不在同一區(qū)域,檢測(cè)新的點(diǎn)是否已經(jīng)被占用if xy_end(pos_id_now)==1%如果該點(diǎn)被占用,說明出現(xiàn)與其它流線太近的情況,則直接停止j=j-1;breakelse%如果沒被占用,則把新點(diǎn)添加上xy_end(pos_id_now)=1;pos_id_last=pos_id_now;endend%順便標(biāo)記xy_startpos_id_s=axis2id(sl_i(j,1),sl_i(j,2),dstart);xy_start(pos_id_s)=1; end sl_i(j:end,:)=[]; endfunction pos_id=axis2id(x,y,distance) N=ceil((0.5-distance/2)/distance)*2+1;%分割的數(shù)量 min_distance=(1-(N-2)*distance)/2;%兩端最小的距離 %x的位置 if x<=min_distancepos_id_x=1; elseif x>=1-min_distancepos_id_x=N; elsepos_id_x=ceil((x-min_distance)/distance)+1; end %y的位置 if y<=min_distancepos_id_y=1; elseif y>=1-min_distancepos_id_y=N; elsepos_id_y=ceil((y-min_distance)/distance)+1; end %xy轉(zhuǎn)ind pos_id=sub2ind([N,N],pos_id_y,pos_id_x); endfunction plot_arrow(xy_lim,xy_ratio,xy_arrow,arrow_color,arrow_width,arrow_direction) %初始化箭頭形狀(歸一化的形狀) arrow_0=[0,0;-1,0.5;-1,-0.5]; %對(duì)方向進(jìn)行歸一化 a_dn=arrow_direction(:)./xy_ratio(:); a_dn=a_dn/sqrt(sum(a_dn.^2)); d=(xy_lim(4)-xy_lim(3)+xy_lim(2)-xy_lim(1))/2; %箭頭對(duì)窗口縮放 arrow_1=arrow_0*arrow_width*0.03*d; %箭頭旋轉(zhuǎn) arrow_2=arrow_1*[a_dn(1),a_dn(2);-a_dn(2),a_dn(1)]; %箭頭變形 xy_ratio_n=xy_ratio/sqrt(sum(xy_ratio.^2));%對(duì)比例尺歸一化 arrow_3=arrow_2.*xy_ratio_n+xy_arrow; fill(arrow_3(:,1),arrow_3(:,2),arrow_color,'EdgeColor','none') end

生成的流線圖為:

2 繪制彩色的短線圖

這個(gè)繪制思路很簡(jiǎn)單,就是取每一個(gè)點(diǎn),然后繪制短流線即可。
代碼如下:

clear clc close all %繪制短線彩色圖(長(zhǎng)短相同)%載入數(shù)據(jù) load wind N=5; xx=x(:,:,N); yy=y(:,:,N); uu=u(:,:,N); vv=v(:,:,N); %xmin=min(xx,[],'all');xmax=max(xx,[],'all'); %ymin=min(yy,[],'all');ymax=max(yy,[],'all'); xmin=min(min(min(x)));xmax=max(max(max(x))); ymin=min(min(min(y)));ymax=max(max(max(y)));V2=sqrt(uu.^2+vv.^2);%計(jì)算速度%繪制變量顏色條 N_color=32; P=V2;%把速度作為變量 %P_max=max(P,[],'all');P_min=min(P,[],'all'); P_max=max(max(max(P)));P_min=min(min(min(P)));Num=numel(xx);%流線數(shù)量for k=1:Numstartx=xx(k);starty=yy(k);sl_i=stream2(xx,yy,uu,vv,startx,starty,[0.1, 20]); sl_sum{k}=sl_i{1};%保存流線 endmcp=colormap(parula(N_color));%生成顏色 [~,~,P_color] = histcounts(P(:),linspace(P_min,P_max,N_color));%顏色索引figure(1) hold on xlim([xmin,xmax]) ylim([ymin,ymax]) xy_ratio = get(gca, 'DataAspectRatio'); xy_ratio = xy_ratio(1:2); xy_lim = axis;for k=1:Numplot(sl_sum{k}(:,1),sl_sum{k}(:,2),'color',mcp(P_color(k),:))if size(sl_sum{k},1)<=1continueend%繪制箭頭arrow_direction=sl_sum{k}(end,:)-sl_sum{k}(end-1,:);plot_arrow(xy_lim,xy_ratio,sl_sum{k}(end,:),mcp(P_color(k),:),0.5,arrow_direction) end hold off caxis([P_min,P_max])%顏色條范圍 colorbarfunction plot_arrow(xy_lim,xy_ratio,xy_arrow,arrow_color,arrow_width,arrow_direction) %初始化箭頭形狀(歸一化的形狀) arrow_0=[0,0;-1,0.5;-1,-0.5]; %對(duì)方向進(jìn)行歸一化 a_dn=arrow_direction(:)./xy_ratio(:); a_dn=a_dn/sqrt(sum(a_dn.^2)); d=(xy_lim(4)-xy_lim(3)+xy_lim(2)-xy_lim(1))/2; %箭頭對(duì)窗口縮放 arrow_1=arrow_0*arrow_width*0.03*d; %箭頭旋轉(zhuǎn) arrow_2=arrow_1*[a_dn(1),a_dn(2);-a_dn(2),a_dn(1)]; %箭頭變形 xy_ratio_n=xy_ratio/sqrt(sum(xy_ratio.^2));%對(duì)比例尺歸一化 arrow_3=arrow_2.*xy_ratio_n+xy_arrow; fill(arrow_3(:,1),arrow_3(:,2),arrow_color,'EdgeColor','none') end

繪制出的效果如下:

這些流線基本都是長(zhǎng)短相同的,因?yàn)閟tream2函數(shù)生成的流線,流線上的點(diǎn)基本是等間距的。如果想做出速度V大的流線箭頭更長(zhǎng),速度小的流線更短的效果,這里用自己編寫的二階RK流線生成算法,進(jìn)行流線的繪制。由于interp2速度太慢,這里用griddedInterpolant函數(shù)進(jìn)行替換,優(yōu)化速度將近10倍(我也沒搞清楚為什么interp2速度這么慢)。

matlab代碼如下:

clear clc close all%繪制彩色短線圖(長(zhǎng)短不同) %載入數(shù)據(jù) load wind N=5; xx=x(:,:,N); yy=y(:,:,N); uu=u(:,:,N); vv=v(:,:,N); %xmin=min(xx,[],'all');xmax=max(xx,[],'all'); %ymin=min(yy,[],'all');ymax=max(yy,[],'all'); xmin=min(min(min(x)));xmax=max(max(max(x))); ymin=min(min(min(y)));ymax=max(max(max(y)));V2=sqrt(uu.^2+vv.^2);%計(jì)算速度 %繪制變量顏色條 N_color=32; P=V2;%指定變量為V2 %P_max=max(P,[],'all');P_min=min(P,[],'all'); P_max=max(max(max(P)));P_min=min(min(min(P)));Num = numel(xx);num_streamline=20;%每條流線上的點(diǎn)數(shù)量 V2_max=max(max(max(V2))); dt=3.0*min([xx(1,2)-xx(1,1),yy(2,1)-yy(1,1)])/V2_max/num_streamline;for k=1:Numstartx=xx(k);starty=yy(k);sl_i=stream2_RK2(xx,yy,uu,vv,startx,starty,dt,num_streamline); sl_sum{k}=sl_i; endmcp=colormap(parula(N_color)); [~,~,P_color] = histcounts(P(:),linspace(P_min,P_max,N_color)); mlw=linspace(0.5,2,N_color);figure(1) hold on xlim([xmin,xmax]) ylim([ymin,ymax]) xy_ratio = get(gca, 'DataAspectRatio'); xy_ratio = xy_ratio(1:2); xy_lim = axis;for k=1:Numplot(sl_sum{k}(:,1),sl_sum{k}(:,2),'color',mcp(P_color(k),:),'linewidth',mlw(P_color(k)))if size(sl_sum{k},1)<=1continueend%繪制箭頭arrow_direction=sl_sum{k}(end,:)-sl_sum{k}(end-1,:);plot_arrow(xy_lim,xy_ratio,sl_sum{k}(end,:),mcp(P_color(k),:),0.5*mlw(P_color(k)),arrow_direction) end hold off caxis([P_min,P_max]) colorbarfunction streamline_i=stream2_RK2(x,y,u,v,startx,starty,dt,N) streamline_i=zeros(N,2); streamline_i(1,:)=[startx,starty]; x_old=startx;y_old=starty; F_u = griddedInterpolant(x',y',u','linear'); F_v = griddedInterpolant(x',y',v','linear'); for k=2:N%利用改進(jìn)歐拉法(或者叫2階Runge-Kutta,預(yù)估校正)%interp2太慢,放棄 % u_K1=interp2(x,y,u,x_old,y_old,'linear')*dt; % v_K1=interp2(x,y,v,x_old,y_old,'linear')*dt;u_K1 = F_u(x_old,y_old)*dt;v_K1 = F_v(x_old,y_old)*dt; % u_K2=interp2(x,y,u,x_old+0.5*u_K1,y_old+0.5*v_K1,'linear')*dt; % v_K2=interp2(x,y,v,x_old+0.5*u_K1,y_old+0.5*v_K1,'linear')*dt;u_K2 = F_u(x_old+0.5*u_K1,y_old+0.5*v_K1)*dt;v_K2 = F_v(x_old+0.5*u_K1,y_old+0.5*v_K1)*dt;x_new=x_old+0.5*(u_K1+u_K2);y_new=y_old+0.5*(v_K1+v_K2);%保存streamline_i(k,:)=[x_new,y_new];x_old=x_new;y_old=y_new;if isnan(x_new) || isnan(y_new)streamline_i(k+1:end,:)=[];breakend end endfunction plot_arrow(xy_lim,xy_ratio,xy_arrow,arrow_color,arrow_width,arrow_direction) %初始化箭頭形狀(歸一化的形狀) arrow_0=[0,0;-1,0.5;-1,-0.5]; %對(duì)方向進(jìn)行歸一化 a_dn=arrow_direction(:)./xy_ratio(:); a_dn=a_dn/sqrt(sum(a_dn.^2)); d=(xy_lim(4)-xy_lim(3)+xy_lim(2)-xy_lim(1))/2; %箭頭對(duì)窗口縮放 arrow_1=arrow_0*arrow_width*0.03*d; %箭頭旋轉(zhuǎn) arrow_2=arrow_1*[a_dn(1),a_dn(2);-a_dn(2),a_dn(1)]; %箭頭變形 xy_ratio_n=xy_ratio/sqrt(sum(xy_ratio.^2));%對(duì)比例尺歸一化 arrow_3=arrow_2.*xy_ratio_n+xy_arrow; fill(arrow_3(:,1),arrow_3(:,2),arrow_color,'EdgeColor','none') end

生成的流線效果如下:

3 繪制彩色的均勻流線

流線的繪制采用第1小節(jié)方法,但是采用控制stream2函數(shù)的點(diǎn)數(shù)來控制流線的長(zhǎng)度(沒有采用第2小節(jié)的RK方法)。顏色的繪制采用第2小節(jié)的方法。

代碼如下:

clear clc close all%繪制彩色長(zhǎng)線圖(長(zhǎng)短不同) %載入數(shù)據(jù) load wind N=5; xx=x(:,:,N); yy=y(:,:,N); uu=u(:,:,N); vv=v(:,:,N); %xmin=min(xx,[],'all');xmax=max(xx,[],'all'); %ymin=min(yy,[],'all');ymax=max(yy,[],'all'); xmin=min(min(min(x)));xmax=max(max(max(x))); ymin=min(min(min(y)));ymax=max(max(max(y)));%流線求解 [streamline_sum,streamline_seed]=my_streamline_mutli(xx,yy,uu,vv,0.03,32);%繪制變量顏色條 V2=sqrt(uu.^2+vv.^2);N_color=16;%流線的顏色種類 P=V2;%把速度作為顏色的變量 %P_max=max(P,[],'all');P_min=min(P,[],'all'); P_max=max(max(max(P)));P_min=min(min(min(P)));mcp=colormap(jet(N_color)); P_seed=interp2(xx,yy,P,streamline_seed(:,1),streamline_seed(:,2));[~,~,P_color] = histcounts(P_seed(:),linspace(P_min,P_max,N_color)); mlw=linspace(0.5,2,N_color);%繪制流線 figure(1) hold on for j=1:size(streamline_sum,2)plot(streamline_sum{j}(:,1),streamline_sum{j}(:,2),'color',mcp(P_color(j),:)...,'linewidth',mlw(P_color(j)));%繪制流線 end xlim([xmin,xmax]) ylim([ymin,ymax]) caxis([P_min,P_max]) colorbar %繪制箭頭 xy_ratio = get(gca, 'DataAspectRatio'); xy_ratio = xy_ratio(1:2); xy_lim = axis; for k=1:size(streamline_sum,2)%繪制箭頭arrow_direction=streamline_sum{k}(end,:)-streamline_sum{k}(end-1,:);plot_arrow(xy_lim,xy_ratio,streamline_sum{k}(end,:),mcp(P_color(k),:),mlw(P_color(k)),arrow_direction) end hold offfunction [streamline_sum,streamline_seed]=my_streamline_mutli(x,y,u,v,dstart,num) %0處理前設(shè)置 %設(shè)置網(wǎng)格密度(01區(qū)間內(nèi)歸一化的長(zhǎng)度) %dstart=0.05;默認(rèn)0.05 dend=0.5*dstart;%xmin=min(x,[],'all');xmax=max(x,[],'all'); %ymin=min(y,[],'all');ymax=max(y,[],'all'); xmin=min(min(min(x)));xmax=max(max(max(x))); ymin=min(min(min(y)));ymax=max(max(max(y)));%歸一化,將流場(chǎng)縮放為0-1區(qū)間的矩形 xn=(x-xmin)/(xmax-xmin); yn=(y-ymin)/(ymax-ymin); un=u/(xmax-xmin); vn=v/(ymax-ymin);num_start=ceil((0.5-dstart/2)/dstart)*2+1; num_end=ceil((0.5-dend/2)/dend)*2+1;%初始化所有網(wǎng)格點(diǎn),0代表可以放置新點(diǎn),1代表已經(jīng)存在原有的點(diǎn) xy_start=zeros(num_start,num_start); xy_end=zeros(num_end,num_end);%將流線劃分為num種,速度越大的流線越長(zhǎng) length_sl=linspace(5,40,num); V2=sqrt(un.^2+vn.^2); V2_max=max(max(max(V2))); V2_min=min(min(min(V2))); %V2_min=min(V2,[],'all');V2_max=max(V2,[],'all');V2_space=linspace(V2_min,V2_max,num+1);%1當(dāng)xy_start內(nèi)還有可放置的新點(diǎn)的位置時(shí),進(jìn)行循環(huán) k=0;%循環(huán)次數(shù),也是流線個(gè)數(shù) while ~all(xy_start,'all')k=k+1;%2隨機(jī)一個(gè)start內(nèi)網(wǎng)格點(diǎn)作為種子點(diǎn)[start_id_y,start_id_x]=find(xy_start==0);randnum=randi(size(start_id_y,1));x_pos_i=id2axis(dstart,start_id_x(randnum,1));y_pos_i=id2axis(dstart,start_id_y(randnum,1));streamline_seed(k,:)=[x_pos_i,y_pos_i];%保存種子點(diǎn)V2_seed=interp2(xn,yn,V2,x_pos_i,y_pos_i);%計(jì)算種子點(diǎn)處的速度[~,~,sl_N] = histcounts(V2_seed,V2_space);num_streamline=round(length_sl(sl_N));%3繪制流線streamline_i_1 = stream2(xn,yn, un, vn,x_pos_i,y_pos_i,[0.1,num_streamline]);streamline_i_2 = stream2(xn,yn,-un,-vn,x_pos_i,y_pos_i,[0.1,num_streamline]);%4以xy_end為標(biāo)準(zhǔn),刪除自相交或間隔太近的點(diǎn)。并順便標(biāo)記xy_end[streamline_i_1,xy_end,xy_start]=delete_self(streamline_i_1{1},xy_end,dend,xy_start,dstart);[streamline_i_2,xy_end,xy_start]=delete_self(streamline_i_2{1},xy_end,dend,xy_start,dstart);%5保存streamline_k=[flipud(streamline_i_2);streamline_i_1(2:end,:)];%新的流線streamline_sum{k}=[xmin+streamline_k(:,1)*(xmax-xmin),ymin+streamline_k(:,2)*(ymax-ymin)];%從歸一化還原 end streamline_seed=[streamline_seed(:,1)*(xmax-xmin)+xmin,streamline_seed(:,2)*(ymax-ymin)+ymin]; endfunction xpoint=id2axis(distance,id) %取網(wǎng)格的中點(diǎn)N=ceil((0.5-distance/2)/distance)*2+1;%分割的數(shù)量 min_distance=(1-(N-2)*distance)/2;%兩端最小的距離 if id==1xpoint=min_distance/2; elseif id==Nxpoint=1-min_distance/2; elsexpoint=min_distance+(id-1.5)*distance; end endfunction [sl_i,xy_end,xy_start]=delete_self(sl_i,xy_end,dend,xy_start,dstart) %sl_i streamline流線,兩列N行形式 N=size(sl_i,1); pos_id_last=axis2id(sl_i(1,1),sl_i(1,2),dend); xy_end(pos_id_last)=1;%第一個(gè)點(diǎn)標(biāo)記%順便標(biāo)記xy_start pos_id_s=axis2id(sl_i(1,1),sl_i(1,2),dstart); xy_start(pos_id_s)=1; for j=2:Npos_id_now=axis2id(sl_i(j,1),sl_i(j,2),dend);if pos_id_now~=pos_id_last%如果現(xiàn)在的點(diǎn)和原有的點(diǎn)在同一區(qū)域,則不管它%如果不在同一區(qū)域,檢測(cè)新的點(diǎn)是否已經(jīng)被占用if xy_end(pos_id_now)==1%如果該點(diǎn)被占用,說明出現(xiàn)與其它流線太近的情況,則直接停止j=j-1;breakelse%如果沒被占用,則把新點(diǎn)添加上xy_end(pos_id_now)=1;pos_id_last=pos_id_now;endend%順便標(biāo)記xy_startpos_id_s=axis2id(sl_i(j,1),sl_i(j,2),dstart);xy_start(pos_id_s)=1; end sl_i(j:end,:)=[]; endfunction pos_id=axis2id(x,y,distance) N=ceil((0.5-distance/2)/distance)*2+1;%分割的數(shù)量 min_distance=(1-(N-2)*distance)/2;%兩端最小的距離%x的位置 if x<=min_distancepos_id_x=1; elseif x>=1-min_distancepos_id_x=N; elsepos_id_x=ceil((x-min_distance)/distance)+1; end%y的位置 if y<=min_distancepos_id_y=1; elseif y>=1-min_distancepos_id_y=N; elsepos_id_y=ceil((y-min_distance)/distance)+1; end%xy轉(zhuǎn)ind pos_id=sub2ind([N,N],pos_id_y,pos_id_x); endfunction plot_arrow(xy_lim,xy_ratio,xy_arrow,arrow_color,arrow_width,arrow_direction) %初始化箭頭形狀(歸一化的形狀) arrow_0=[0,0;-1,0.5;-1,-0.5]; %對(duì)方向進(jìn)行歸一化 a_dn=arrow_direction(:)./xy_ratio(:); a_dn=a_dn/sqrt(sum(a_dn.^2)); d=(xy_lim(4)-xy_lim(3)+xy_lim(2)-xy_lim(1))/2; %箭頭對(duì)窗口縮放 arrow_1=arrow_0*arrow_width*0.03*d; %箭頭旋轉(zhuǎn) arrow_2=arrow_1*[a_dn(1),a_dn(2);-a_dn(2),a_dn(1)]; %箭頭變形 xy_ratio_n=xy_ratio/sqrt(sum(xy_ratio.^2));%對(duì)比例尺歸一化 %arrow_3=arrow_2.*xy_ratio_n+xy_arrow;%由于兼容性問題,更改為下面語句 xy_ratio_n2=ones(3,1)*xy_ratio_n; arrow_3=arrow_2.*xy_ratio_n2+xy_arrow; fill(arrow_3(:,1),arrow_3(:,2),arrow_color,'EdgeColor','none') end

4 運(yùn)動(dòng)的彩色箭頭流線圖

一開始我本來是打算用第2小節(jié)的方式繪制一個(gè)動(dòng)圖,但是發(fā)現(xiàn)原來的箭頭(粒子)離開之后,沒有新的箭頭進(jìn)行補(bǔ)充。如果明確流場(chǎng)的入口和出口,可以采用出口離開多少粒子,就在入口添加多少粒子的方式,維持流場(chǎng)內(nèi)數(shù)量恒定。但是對(duì)于復(fù)雜流場(chǎng)則不容易做到這點(diǎn)。

因此,我采用了均勻流線的方式,每次運(yùn)動(dòng)完成后,刪除離得太近的流線,然后在空白處補(bǔ)充新的流線。效果如下:

這樣做的優(yōu)點(diǎn)在于每一幀都滿足均勻流線的要求,而且保證絕大多數(shù)流線不會(huì)憑空的產(chǎn)生和消失,即使是消失也是逐漸消失。

代碼如下:

clear clc close all %繪制彩色長(zhǎng)線圖(長(zhǎng)短不同,而且是會(huì)動(dòng)的那種gif圖)%載入數(shù)據(jù) load wind N=5; xx=x(:,:,N); yy=y(:,:,N); uu=u(:,:,N); vv=v(:,:,N); %xmin=min(xx,[],'all');xmax=max(xx,[],'all'); %ymin=min(yy,[],'all');ymax=max(yy,[],'all'); xmin=min(min(min(x)));xmax=max(max(max(x))); ymin=min(min(min(y)));ymax=max(max(max(y)));dstart=0.03;%網(wǎng)格寬度%流線求解 [streamline_sum,streamline_seed]=my_streamline_mutli(xx,yy,uu,vv,dstart,32);%繪制變量顏色條 V2=sqrt(uu.^2+vv.^2);N_color=16;%流線的顏色種類 P=V2;%把速度作為顏色的變量 %P_max=max(P,[],'all');P_min=min(P,[],'all'); P_max=max(max(max(P)));P_min=min(min(min(P)));mcp=colormap(jet(N_color)); P_seed=interp2(xx,yy,P,streamline_seed(:,1),streamline_seed(:,2));[~,~,P_color] = histcounts(P_seed(:),linspace(P_min,P_max,N_color)); mlw=linspace(0.5,2,N_color);%繪制流線 figure(1) hold on for j=1:size(streamline_sum,2)plot(streamline_sum{j}(:,1),streamline_sum{j}(:,2),'color',mcp(P_color(j),:)...,'linewidth',mlw(P_color(j)));%繪制流線 end xlim([xmin,xmax]) ylim([ymin,ymax])%固定圖窗大小 caxis([P_min,P_max])%設(shè)定顏色條范圍 colorbar%顯示色條xy_ratio = get(gca, 'DataAspectRatio'); xy_ratio = xy_ratio(1:2); xy_lim = axis; for k=1:size(streamline_sum,2)%繪制箭頭if size(streamline_sum{k},1)<=2continueendarrow_direction=streamline_sum{k}(end,:)-streamline_sum{k}(end-1,:);plot_arrow(xy_lim,xy_ratio,streamline_sum{k}(end,:),mcp(P_color(k),:),mlw(P_color(k)),arrow_direction) end hold off%保存gif圖 frame=getframe(gca); im=frame2im(frame); [I,map]=rgb2ind(im,36); imwrite(I,map,'temp.gif','gif', 'Loopcount',inf,'DelayTime',0.08);%第一張%后續(xù)的變化 for k=1:60%延長(zhǎng)dt秒之后的圖像streamline_sum=my_streamline_time(xx,yy,uu,vv,streamline_sum,dstart,0.05);figure(1)clfhold onfor j=1:size(streamline_sum,2)sl_i=streamline_sum{j};if isempty(sl_i)continueendL_sl_i=size(sl_i,1);P_sl=interp2(xx,yy,P,sl_i(round(L_sl_i/2),1),sl_i(round(L_sl_i/2),2));[~,~,P_color] = histcounts(P_sl,linspace(P_min,P_max,N_color));P_color(P_color==0)=1;plot(sl_i(:,1),sl_i(:,2),'color',mcp(P_color,:)...,'linewidth',mlw(P_color));%繪制流線%繪制箭頭if size(sl_i,1)<=2continueendarrow_direction=sl_i(end,:)-sl_i(end-1,:);plot_arrow(xy_lim,xy_ratio,sl_i(end,:),mcp(P_color,:),mlw(P_color),arrow_direction)endhold offxlim([xmin,xmax])ylim([ymin,ymax])%固定圖窗大小caxis([P_min,P_max])%設(shè)定顏色條范圍colorbar%顯示色條pause(0.5)%保存gif圖frame=getframe(gca);im=frame2im(frame);[I,map]=rgb2ind(im,36);imwrite(I,map,'temp.gif','gif','WriteMode','append','DelayTime',0.08); endfunction [streamline_sum,streamline_seed]=my_streamline_mutli(x,y,u,v,dstart,num) %繪制流線 %0處理前設(shè)置%設(shè)置網(wǎng)格密度(01區(qū)間內(nèi)歸一化的長(zhǎng)度) %dstart=0.05;默認(rèn)0.05 dend=0.5*dstart; %xmin=min(x,[],'all');xmax=max(x,[],'all'); %ymin=min(y,[],'all');ymax=max(y,[],'all'); xmin=min(min(min(x)));xmax=max(max(max(x))); ymin=min(min(min(y)));ymax=max(max(max(y)));%歸一化,將流場(chǎng)縮放為0-1區(qū)間的矩形 xn=(x-xmin)/(xmax-xmin); yn=(y-ymin)/(ymax-ymin); un=u/(xmax-xmin); vn=v/(ymax-ymin); F_u = griddedInterpolant(xn',yn',un','linear'); F_v = griddedInterpolant(xn',yn',vn','linear'); F_u_n = griddedInterpolant(xn',yn',-un','linear'); F_v_n = griddedInterpolant(xn',yn',-vn','linear'); num_start=ceil((0.5-dstart/2)/dstart)*2+1; num_end=ceil((0.5-dend/2)/dend)*2+1;%初始化所有網(wǎng)格點(diǎn),0代表可以放置新點(diǎn),1代表已經(jīng)存在原有的點(diǎn) xy_start=zeros(num_start,num_start); xy_end=zeros(num_end,num_end); %將流線劃分為num種,速度越大的流線越長(zhǎng) length_sl=linspace(5,40,num); V2=sqrt(un.^2+vn.^2); %V2_min=min(V2,[],'all');V2_max=max(V2,[],'all'); V2_max=max(max(max(V2))); V2_min=min(min(min(V2)));V2_space=linspace(V2_min,V2_max,num+1);%1當(dāng)xy_start內(nèi)還有可放置的新點(diǎn)的位置時(shí),進(jìn)行循環(huán) k=0;%循環(huán)次數(shù),也是流線個(gè)數(shù) while ~all(xy_start,'all')k=k+1;%2隨機(jī)一個(gè)start內(nèi)網(wǎng)格點(diǎn)作為種子點(diǎn)[start_id_y,start_id_x]=find(xy_start==0);randnum=randi(size(start_id_y,1));x_pos_i=id2axis(dstart,start_id_x(randnum,1));y_pos_i=id2axis(dstart,start_id_y(randnum,1));streamline_seed(k,:)=[x_pos_i,y_pos_i];%保存種子點(diǎn)%3繪制流線streamline_i_1=stream2_RK2(F_u ,F_v ,x_pos_i,y_pos_i,0.01,20);streamline_i_2=stream2_RK2(F_u_n,F_v_n,x_pos_i,y_pos_i,0.01,20);%4以xy_end為標(biāo)準(zhǔn),刪除自相交或間隔太近的點(diǎn)。并順便標(biāo)記xy_end[streamline_i_1,xy_end,xy_start]=delete_self(streamline_i_1,xy_end,dend,xy_start,dstart);[streamline_i_2,xy_end,xy_start]=delete_self(streamline_i_2,xy_end,dend,xy_start,dstart);%5保存streamline_k=[flipud(streamline_i_2);streamline_i_1(2:end,:)];%新的流線streamline_sum{k}=[xmin+streamline_k(:,1)*(xmax-xmin),ymin+streamline_k(:,2)*(ymax-ymin)];%從歸一化還原 end streamline_seed=[streamline_seed(:,1)*(xmax-xmin)+xmin,streamline_seed(:,2)*(ymax-ymin)+ymin]; endfunction xpoint=id2axis(distance,id) %取網(wǎng)格的中點(diǎn) N=ceil((0.5-distance/2)/distance)*2+1;%分割的數(shù)量 min_distance=(1-(N-2)*distance)/2;%兩端最小的距離 if id==1xpoint=min_distance/2; elseif id==Nxpoint=1-min_distance/2; elsexpoint=min_distance+(id-1.5)*distance; end endfunction [sl_i,xy_end,xy_start]=delete_self(sl_i,xy_end,dend,xy_start,dstart) %sl_i streamline流線,兩列N行形式 N=size(sl_i,1); pos_id_last=axis2id(sl_i(1,1),sl_i(1,2),dend); xy_end(pos_id_last)=1;%第一個(gè)點(diǎn)標(biāo)記%順便標(biāo)記xy_start pos_id_s=axis2id(sl_i(1,1),sl_i(1,2),dstart); xy_start(pos_id_s)=1; for j=2:Npos_id_now=axis2id(sl_i(j,1),sl_i(j,2),dend);if pos_id_now~=pos_id_last%如果現(xiàn)在的點(diǎn)和原有的點(diǎn)在同一區(qū)域,則不管它%如果不在同一區(qū)域,檢測(cè)新的點(diǎn)是否已經(jīng)被占用if xy_end(pos_id_now)==1%如果該點(diǎn)被占用,說明出現(xiàn)與其它流線太近的情況,則直接停止j=j-1;breakelse%如果沒被占用,則把新點(diǎn)添加上xy_end(pos_id_now)=1;pos_id_last=pos_id_now;endend%順便標(biāo)記xy_startpos_id_s=axis2id(sl_i(j,1),sl_i(j,2),dstart);xy_start(pos_id_s)=1; end sl_i(j:end,:)=[];%刪除停止之后剩余的流線 endfunction pos_id=axis2id(x,y,distance) N=ceil((0.5-distance/2)/distance)*2+1;%分割的數(shù)量 min_distance=(1-(N-2)*distance)/2;%兩端最小的距離 %x的位置 if x<=min_distancepos_id_x=1; elseif x>=1-min_distancepos_id_x=N; elsepos_id_x=ceil((x-min_distance)/distance)+1; end %y的位置 if y<=min_distancepos_id_y=1; elseif y>=1-min_distancepos_id_y=N; elsepos_id_y=ceil((y-min_distance)/distance)+1; end %xy轉(zhuǎn)ind pos_id=sub2ind([N,N],pos_id_y,pos_id_x); endfunction plot_arrow(xy_lim,xy_ratio,xy_arrow,arrow_color,arrow_width,arrow_direction) %初始化箭頭形狀(歸一化的形狀) arrow_0=[0,0;-1,0.5;-1,-0.5]; %對(duì)方向進(jìn)行歸一化 a_dn=arrow_direction(:)./xy_ratio(:); a_dn=a_dn/sqrt(sum(a_dn.^2)); d=(xy_lim(4)-xy_lim(3)+xy_lim(2)-xy_lim(1))/2; %箭頭對(duì)窗口縮放 arrow_1=arrow_0*arrow_width*0.03*d; %箭頭旋轉(zhuǎn) arrow_2=arrow_1*[a_dn(1),a_dn(2);-a_dn(2),a_dn(1)]; %箭頭變形 xy_ratio_n=xy_ratio/sqrt(sum(xy_ratio.^2));%對(duì)比例尺歸一化 %arrow_3=arrow_2.*xy_ratio_n+xy_arrow;%由于兼容性問題,更改為下面語句 xy_ratio_n2=ones(3,1)*xy_ratio_n; arrow_3=arrow_2.*xy_ratio_n2+xy_arrow; fill(arrow_3(:,1),arrow_3(:,2),arrow_color,'EdgeColor','none') endfunction streamline_i=stream2_RK2(F_u,F_v,startx,starty,dt,N) %繪制流線,采用RK方法 %利用改進(jìn)歐拉法(或者叫2階Runge-Kutta,預(yù)估校正) streamline_i=zeros(N,2); streamline_i(1,:)=[startx,starty]; x_old=startx;y_old=starty; for k=2:Nu_K1 = F_u(x_old,y_old)*dt;v_K1 = F_v(x_old,y_old)*dt;u_K2 = F_u(x_old+0.5*u_K1,y_old+0.5*v_K1)*dt;v_K2 = F_v(x_old+0.5*u_K1,y_old+0.5*v_K1)*dt;x_new=x_old+0.5*(u_K1+u_K2);y_new=y_old+0.5*(v_K1+v_K2);%保存streamline_i(k,:)=[x_new,y_new];x_old=x_new;y_old=y_new;if isnan(x_new) || isnan(y_new)streamline_i(k+1:end,:)=[];breakendif 0>x_new || 0>y_new || x_new>1 || y_new>1%由于插值結(jié)果都是在0-1區(qū)間內(nèi)的歸一化數(shù)據(jù),所以超出邊界的刪除。streamline_i(k+1:end,:)=[];breakend end endfunction streamline_sum=my_streamline_time(x,y,u,v,streamline_sum,dstart,dt) %0處理前設(shè)置 %xmin=min(x,[],'all');xmax=max(x,[],'all'); %ymin=min(y,[],'all');ymax=max(y,[],'all'); xmin=min(min(min(x)));xmax=max(max(max(x))); ymin=min(min(min(y)));ymax=max(max(max(y)));%歸一化,將流場(chǎng)縮放為0-1區(qū)間的矩形 xn=(x-xmin)/(xmax-xmin); yn=(y-ymin)/(ymax-ymin); un=u/(xmax-xmin); vn=v/(ymax-ymin);streamline_sum(cellfun(@isempty,streamline_sum))=[];%刪除空數(shù)組 %排序,優(yōu)先保證長(zhǎng)線段不被刪除 L_sl=zeros(length(streamline_sum),1); for k=1:length(streamline_sum)L_sl(k)=sum((streamline_sum{k}(end,:)-streamline_sum{k}(1,:)).^2); end [~,I]=sort(L_sl,'descend' ); streamline_sum(:,[1:length(streamline_sum)])=streamline_sum(:,[I]);%設(shè)置網(wǎng)格密度(01區(qū)間內(nèi)歸一化的長(zhǎng)度) %dstart=0.05;默認(rèn)0.05 dend=0.5*dstart;%每一條流線隨時(shí)間的變化 F_u = griddedInterpolant(xn',yn',un','linear'); F_v = griddedInterpolant(xn',yn',vn','linear'); F_u_n = griddedInterpolant(xn',yn',-un','linear'); F_v_n = griddedInterpolant(xn',yn',-vn','linear');num_start=ceil((0.5-dstart/2)/dstart)*2+1; num_end=ceil((0.5-dend/2)/dend)*2+1;for k=1:length(streamline_sum)sl_i=streamline_sum{k};if isempty(sl_i)continueendif any(isnan(sl_i(end,:)))sl_i(1,:)=[];elsesl_i(1,:)=[];if isempty(sl_i)continueendsl_i_n=(sl_i-[xmin,ymin])./[(xmax-xmin),(ymax-ymin)];%歸一化流線sl_i_now=stream2_RK2(F_u ,F_v ,sl_i_n(end,1),sl_i_n(end,2),dt,2);sl_i_n=[sl_i_n;sl_i_now(2,:)];sl_i=sl_i_n.*[(xmax-xmin),(ymax-ymin)]+[xmin,ymin];%反歸一化endstreamline_sum{k}=sl_i; end streamline_sum(cellfun(@isempty,streamline_sum))=[];%刪除空數(shù)組xy_start=zeros(num_start,num_start); xy_end=zeros(num_end,num_end); for k=1:length(streamline_sum)sl_i=streamline_sum{k};sl_i=flipud(sl_i);%刪除尾部,保留頭部sl_i_n=(sl_i-[xmin,ymin])./[(xmax-xmin),(ymax-ymin)];%歸一化流線%以xy_end為標(biāo)準(zhǔn),刪除自相交或間隔太近的點(diǎn)。并順便標(biāo)記xy_end[sl_i_n,xy_end,xy_start]=delete_self(sl_i_n,xy_end,dend,xy_start,dstart);sl_i=sl_i_n.*[(xmax-xmin),(ymax-ymin)]+[xmin,ymin];%反歸一化sl_i=flipud(sl_i);streamline_sum{k}=sl_i; end %之后還有一些區(qū)域太空,所以重新在那些區(qū)域生成流線 while ~all(xy_start,'all')k=k+1;%2隨機(jī)一個(gè)start內(nèi)網(wǎng)格點(diǎn)作為種子點(diǎn)[start_id_y,start_id_x]=find(xy_start==0);randnum=randi(size(start_id_y,1));x_pos_i=id2axis(dstart,start_id_x(randnum,1));y_pos_i=id2axis(dstart,start_id_y(randnum,1));%3繪制流線streamline_i_1=stream2_RK2(F_u ,F_v ,x_pos_i,y_pos_i,0.01,20);streamline_i_2=stream2_RK2(F_u_n,F_v_n,x_pos_i,y_pos_i,0.01,20);%4以xy_end為標(biāo)準(zhǔn),刪除自相交或間隔太近的點(diǎn)。并順便標(biāo)記xy_end[streamline_i_1,xy_end,xy_start]=delete_self(streamline_i_1,xy_end,dend,xy_start,dstart);[streamline_i_2,xy_end,xy_start]=delete_self(streamline_i_2,xy_end,dend,xy_start,dstart);%5保存sl_i_n=[flipud(streamline_i_2);streamline_i_1(2:end,:)];%新的流線sl_i=sl_i_n.*[(xmax-xmin),(ymax-ymin)]+[xmin,ymin];%反歸一化streamline_sum{k}=sl_i; endend

總結(jié)

以上是生活随笔為你收集整理的利用matlab绘制二维均匀流线和向量场的全部?jī)?nèi)容,希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯(cuò),歡迎將生活随笔推薦給好友。

国产精品网站 | 国产精品久久久区三区天天噜 | 成人免费影院 | 国产专区一 | 亚洲一二三久久 | 国产精品国产亚洲精品看不卡 | 色姑娘综合 | 日日干网| 免费一级特黄毛大片 | 日日夜精品 | 亚洲黄网址| 曰本三级在线 | 亚洲精品视频国产 | 人人干人人搞 | 色资源网在线观看 | 91精品国产乱码在线观看 | 在线观看av免费观看 | 夜添久久精品亚洲国产精品 | 久久色在线观看 | 在线观看日韩av | 黄网av在线 | 国产精品国内免费一区二区三区 | 欧美日本国产在线观看 | 天天操夜夜逼 | 亚洲三级国产 | 福利av影院 | 婷婷中文字幕 | 麻豆视频国产 | 91私密视频 | 欧美专区亚洲专区 | 96av在线视频 | 亚洲视频免费在线看 | 99成人精品| 超碰官网 | 日本视频不卡 | 99精品国产一区二区三区麻豆 | 丝袜一区在线 | 亚洲九九九在线观看 | 亚洲jizzjizz日本少妇 | 国产伦理一区二区三区 | 日韩欧美成| 亚洲欧美在线综合 | www.国产在线视频 | 六月丁香久久 | 日韩精品中文字幕有码 | 91久久久国产精品 | 久久九九国产视频 | 成人影片在线播放 | 97色婷婷成人综合在线观看 | 99色资源 | 欧美精品乱码久久久久久按摩 | 91日韩国产| 亚洲动漫在线观看 | 免费av在线网站 | 99视频精品免费观看, | 国产成人精品一区一区一区 | 国产美女精品在线 | 国产在线欧美日韩 | 黄色avwww| 在线免费亚洲 | 国产99久久久欧美黑人 | 日韩精品一区二区免费 | 亚洲欧美综合 | 探花视频免费观看高清视频 | 国产高清视频在线播放 | 97视频在线观看视频免费视频 | 韩日色视频 | 日韩视频在线观看视频 | 久久精品这里热有精品 | 国产裸体永久免费视频网站 | 精品亚洲免费视频 | 夜夜骑首页| 日韩一三区 | 综合久久久久久久久 | 国产免费亚洲高清 | 深夜免费小视频 | 一区二区三区电影 | 在线视频在线观看 | 日韩在线视频观看免费 | 国产亚洲欧美精品久久久久久 | 国产福利91精品一区二区三区 | 日韩毛片一区 | 天天射天天色天天干 | 国产视频在线看 | 国产精品国产亚洲精品看不卡 | 久久久av免费 | 中文字幕精| 人人爽人人看 | 中文字幕在线播放一区 | www在线观看视频 | 久草在线视频在线观看 | av在线亚洲天堂 | 欧美成人精品欧美一级乱黄 | 91视频xxxx| www.在线观看视频 | 国产一区二区三区在线 | 91麻豆精品国产自产 | 国产视频丨精品|在线观看 国产精品久久久久久久久久久久午夜 | 99免费在线播放99久久免费 | 久草在线最新 | 久久免费观看少妇a级毛片 久久久久成人免费 | 国产xxxxx在线观看 | 欧美成人精品xxx | 免费激情网 | 深夜福利视频一区二区 | 色资源在线 | 精品99免费视频 | 免费久久片 | 国产美女在线免费观看 | 99产精品成人啪免费网站 | a极黄色片 | 国产精品 欧美 日韩 | 国产精品久久久久久久久免费 | www.天天干.com | 日本系列中文字幕 | 亚洲第一久久久 | 骄小bbw搡bbbb揉bbbb | 中文字幕在线免费97 | aav在线| 天天天色综合a | 国产精品9999久久久久仙踪林 | 久久久久成人精品免费播放动漫 | 国产aaa免费视频 | 欧美激情精品久久久久久免费印度 | 日韩在线视频网址 | 精品一区二区日韩 | 国产精品a久久久久 | 亚洲视频免费在线 | 亚洲欧美日韩一区二区三区在线观看 | 国产91粉嫩白浆在线观看 | 区一区二区三区中文字幕 | 中文字幕av影院 | 最近中文字幕免费观看 | 亚洲乱码精品久久久久 | 中文字幕日韩电影 | 99精品免费网 | 久久久伊人网 | 伊人天天干 | 国产精品久久久久久久久久久免费看 | 99久久久久国产精品免费 | 国产一级在线视频 | 久久国产露脸精品国产 | 在线观看免费国产小视频 | 中文字幕电影高清在线观看 | 在线免费性生活片 | 日本少妇视频 | 美女网站在线观看 | 青草草在线视频 | 中文字幕在线看视频国产中文版 | 中国一级特黄毛片大片久久 | 日本字幕网| 国产一区二区在线观看视频 | 国产打女人屁股调教97 | 一区二区视频免费在线观看 | 午夜精品久久久久久99热明星 | 国产不卡视频在线播放 | 九九视频这里只有精品 | 最近更新好看的中文字幕 | 91成品人影院 | 五月激情婷婷丁香 | 亚洲精品久久久久中文字幕m男 | 久久久久久高潮国产精品视 | 国产一级在线播放 | 欧美精品亚洲精品日韩精品 | 成人动图| 黄视频色网站 | 在线最新av | 在线观看你懂的网址 | 久久人人97超碰精品888 | 欧美中文字幕第一页 | 精品在线视频播放 | 正在播放国产精品 | 97在线观看 | 国产一区在线观看视频 | 久久久精品 一区二区三区 国产99视频在线观看 | 激情五月婷婷激情 | 美女网色| 国产手机视频在线播放 | 69亚洲精品 | 久久www免费人成看片高清 | 成人福利在线播放 | 国产免费午夜 | 欧美性春潮 | 欧美日韩高清在线观看 | 免费亚洲视频在线观看 | 久久综合免费视频 | 97色在线| 午夜精品剧场 | 国产精品久久久久久久久久久久午夜片 | 99精品在线视频观看 | 亚洲日本va在线观看 | 国产一级片免费视频 | 亚洲精品在线观看的 | 日韩簧片在线观看 | 日韩在线观看网站 | 日韩在线播放欧美字幕 | 手机成人av在线 | 在线影院av | 男女视频久久久 | 中文字幕第一页在线视频 | 一区二区三区国 | 99热这里只有精品1 av中文字幕日韩 | 婷婷综合网 | 久久婷婷丁香 | 香蕉网在线 | 99久久精品费精品 | 久久韩国免费视频 | 亚洲三级毛片 | 欧美性生爱| 久久精品成人热国产成 | 日日夜夜中文字幕 | 午夜精品久久久久久99热明星 | 九九久| 99国产一区二区三精品乱码 | 超碰人人射 | 日日躁夜夜躁aaaaxxxx | 亚州av成人| 国产精品露脸在线 | 国产123区在线观看 国产精品麻豆91 | 四虎免费在线观看 | 五月天丁香视频 | 久久你懂得 | 亚洲欧洲一区二区在线观看 | 九色精品在线 | av大全在线观看 | 97超碰.com | 久久精品影片 | 伊人久久在线观看 | 91麻豆精品国产91久久久无限制版 | а天堂中文最新一区二区三区 | 97精品久久 | 日韩在线三级 | 午夜精品久久久久久久99热影院 | 午夜精品在线看 | 亚洲第一久久久 | 美女福利视频网 | 免费欧美高清视频 | 国产日韩精品一区二区三区 | www.av中文字幕.com | 精品一区精品二区 | 九九热只有这里有精品 | 婷婷丁香七月 | 超碰在线人人爱 | 国产香蕉97碰碰久久人人 | 国产尤物一区二区三区 | 午夜精品福利一区二区 | 亚洲精品久久久蜜桃直播 | 2019久久精品 | 91成人在线视频观看 | 丁香六月在线观看 | 2023av| av电影中文| 久久综合99| 免费观看国产精品视频 | 黄色影院在线播放 | 日韩中文字幕免费在线播放 | 黄色av一级片| 色偷偷男人的天堂av | 国产亚洲成av片在线观看 | 激情视频91| 2022久久国产露脸精品国产 | 国产中文字幕大全 | 国产亚洲精品久久久久5区 成人h电影在线观看 | 亚洲国产精品一区二区尤物区 | 天天色天天上天天操 | 国产91学生粉嫩喷水 | 9797在线看片亚洲精品 | 久久免费视频6 | 精品一区二区在线免费观看 | 国产成人精品亚洲日本在线观看 | 国产精品亚州 | 国产在线一线 | 成人久久综合 | 狠狠操电影网 | 成人h动漫精品一区二 | 国产在线精品一区 | 免费看毛片在线 | 日韩午夜电影 | 国产欧美精品一区二区三区 | 亚洲毛片在线观看. | 黄色片网站av | 色狠狠久久av五月综合 | 日韩亚洲国产精品 | 国产精品一区二区白浆 | 久章草在线观看 | 97超碰在线视 | 最新国产精品久久精品 | 国产96精品 | 日本不卡一区二区 | 国产在线观看h | 色多多污污在线观看 | 91视频a | 日韩av影视 | 一区二区三区日韩精品 | 国产一区二区视频在线播放 | 五月天色站 | 在线看国产视频 | 亚洲成成品网站 | 在线a亚洲视频播放在线观看 | 亚洲欧美国产视频 | 91精品久久久久久久久久入口 | 日韩高清一| 天天操天天干天天 | 亚洲成人av片在线观看 | 不卡日韩av| 欧美日韩精品在线视频 | 久久久久久久久久久久国产精品 | 欧美精品乱码久久久久 | 91手机在线看片 | 97在线看片 | 中文字幕 成人 | 国产美女免费视频 | 天天操比 | 久久手机免费观看 | 国产涩图 | 国产精品久久久久久久7电影 | 国产老妇av| 日韩av午夜 | 欧美日韩在线视频免费 | 久久亚洲福利视频 | 永久免费精品视频 | 日韩av手机在线观看 | 精品亚洲男同gayvideo网站 | 91免费版在线 | 精品亚洲va在线va天堂资源站 | 欧美精品v国产精品v日韩精品 | 日韩中文字幕a | 久久在线观看视频 | 欧美影院久久 | 亚洲国产偷 | 久久电影网站中文字幕 | 99视频精品全国免费 | 免费成人在线网站 | 欧美日韩国产成人 | 中文国产在线观看 | 国语精品免费视频 | 国产成人久久精品 | 超级碰碰免费视频 | 丁香婷婷激情 | 九色在线视频 | 久久99亚洲精品久久 | 天天综合天天做天天综合 | 伊人中文在线 | 久久免费成人网 | 亚洲va韩国va欧美va精四季 | 国产一线二线三线性视频 | 欧美精品在线免费 | 三级在线播放视频 | 国内精品久久久久久久久久清纯 | 91热精品视频| 在线观看亚洲专区 | 人人爽夜夜爽 | 国产福利一区二区三区在线观看 | 亚洲精品综合在线 | 免费网站黄 | 久久伊99综合婷婷久久伊 | www.精选视频.com | 美女网站视频久久 | 久草在线高清 | 美女搞黄国产视频网站 | 久久乐九色婷婷综合色狠狠182 | 日本三级中文字幕在线观看 | 99精品国产一区二区三区不卡 | 亚洲激情在线观看 | 色婷婷成人网 | 日韩av图片| 粉嫩av一区二区三区入口 | 黄色毛片在线观看 | 亚洲国产成人av网 | 久久综合9988久久爱 | 午夜美女wwww| 久章草在线观看 | 日韩欧美视频一区二区 | 九九有精品 | 日韩精品久久一区二区 | 亚洲激情六月 | 日韩高清一区在线 | 久久精品一区八戒影视 | 精品一区二区视频 | 欧美激情精品久久久久久免费 | 日韩电影在线观看中文字幕 | 一本一本久久a久久精品综合妖精 | 成人污视频在线观看 | 丁香影院在线 | 日韩欧美高清在线观看 | 婷婷丁香在线观看 | 久草网在线视频 | 久久久影视 | 免费av福利 | 欧美黄色成人 | 久久久久久久久久久网站 | 91色网址 | 日韩一区二区三区观看 | 国产视频在线观看一区 | 婷婷色综合 | 美女在线观看av | 91人人网| 色噜噜在线观看视频 | 亚洲精品一区二区网址 | 99在线视频免费观看 | 99这里只有 | 99精品在线免费 | 天天草天天干天天射 | 97福利在线| 中文字幕亚洲在线观看 | 开心激情久久 | 国产精品v欧美精品v日韩 | 在线观看香蕉视频 | 人人干干人人 | 日韩欧美一区二区三区在线 | 国产精品久久综合 | 国产一级a毛片视频爆浆 | 久久国产精品99久久人人澡 | 夜夜夜影院 | 国产综合小视频 | a国产精品| 999久久久国产精品 高清av免费观看 | 97精品国产97久久久久久粉红 | 99精品热| 人人射人人爽 | 99久久er热在这里只有精品66 | 好看的国产精品视频 | 国产一区二区在线观看视频 | 久久久久亚洲精品男人的天堂 | 91久久丝袜国产露脸动漫 | 少妇啪啪av入口 | 精品久久久久国产免费第一页 | 日韩二区在线观看 | 丁香综合五月 | 天天躁日日| 激情丁香 | 五月天婷婷狠狠 | 国产专区精品 | 2021久久| 国产精久久久久久妇女av | 国产精品久久久久久久久岛 | 日韩xxxxxxxxx | 9ⅰ精品久久久久久久久中文字幕 | 欧美成人亚洲成人 | 亚洲最新av在线网站 | 日韩电影一区二区在线 | 国产免费观看视频 | 88av视频| 亚州精品一二三区 | 99久久日韩精品视频免费在线观看 | 久草精品视频 | 欧美综合色| 亚洲欧美日韩精品一区二区 | 2023国产精品自产拍在线观看 | 久久99视频 | 欧美日韩p片| 色网站视频 | 激情久久综合 | 97精品超碰一区二区三区 | 99产精品成人啪免费网站 | 在线日韩视频 | 女人18片 | 日韩欧美视频在线免费观看 | 麻豆国产网站 | 九九久久久| 在线观看成人国产 | 国产视频资源 | 久久精品美女视频网站 | 国产a网站 | 成人a级网站 | 91精品久久久久久久久 | 久久激情婷婷 | 国产乱对白刺激视频不卡 | 开心激情综合网 | 成人午夜电影久久影院 | 五月婷婷丁香在线观看 | 日韩中出在线 | 色婷婷www | 91av手机在线观看 | 97视频人人澡人人爽 | 狠狠干网址 | 一区二区三区免费播放 | 日韩精品久久久久久中文字幕8 | 最近能播放的中文字幕 | 玖玖爱在线观看 | 国产中文| 亚洲四虎在线 | 一区中文字幕在线观看 | www黄色大片 | 欧洲黄色片 | 91入口在线观看 | 国产视频一二区 | 免费色视频在线 | 国产成人久久精品一区二区三区 | 久久国产免费 | 色偷偷888欧美精品久久久 | 国产色在线视频 | 久草在线视频资源 | 日韩欧美一区二区三区视频 | 国产流白浆高潮在线观看 | 啪一啪在线 | 亚洲免费视频观看 | 日本中文字幕观看 | 日韩高清www| 日韩精品视频在线免费观看 | 中文字幕av专区 | 精品国产乱码久久久久久1区二区 | 在线看国产一区 | 字幕网av | av免费片| 黄色小网站免费看 | 97香蕉超级碰碰久久免费软件 | 91麻豆产精品久久久久久 | 欧美大香线蕉线伊人久久 | 国产成人免费在线 | 在线激情av电影 | www91在线观看 | 中文日韩在线视频 | 午夜精品婷婷 | 久久综合桃花 | 国产精品久久久久三级 | av 一区二区三区 | 国产永久免费高清在线观看视频 | 午夜精品久久久久久中宇69 | 欧美久久久久久久 | 米奇影视7777| 992tv人人草| 日韩h在线观看 | 日韩精品不卡在线 | 992tv成人免费看片 | 日韩av免费一区二区 | 国产看片免费 | www.久久精品视频 | 国产一区二区三区高清播放 | 久久久久久综合 | 亚洲视频精品在线 | 丁香六月激情 | 欧洲视频一区 | 国产专区视频在线 | 久久久久日本精品一区二区三区 | 91视频在线观看免费 | 国产乱视频 | 亚洲国产成人在线 | 一级免费看| 亚洲黄色一级电影 | 亚洲国产精品免费 | 免费电影播放 | 亚洲理论片在线观看 | 亚洲精品福利视频 | 亚洲一级理论片 | 亚洲乱码精品久久久 | 久久美女高清视频 | 亚洲欧洲精品久久 | 狠狠色噜噜狠狠狠狠 | 国产一级免费观看视频 | 日日夜夜人人精品 | 欧美精品在线观看免费 | 欧美一级片在线播放 | 视频国产一区二区三区 | 国内精品久久久久国产 | 国内精品视频在线 | 久久久黄视频 | av电影在线免费观看 | 五月婷av | 特级xxxxx欧美| 国产视频18| 99久久99久久精品国产片果冰 | 草久在线观看 | 国产青青青| 天堂在线一区二区三区 | 亚洲精品xxxx | 伊人电影天堂 | 亚洲欧洲中文日韩久久av乱码 | 久久超碰免费 | 中文字幕av免费在线观看 | 天堂网一区二区 | 毛片精品免费在线观看 | 91免费高清 | 精品一区二区三区久久久 | 国产永久免费高清在线观看视频 | 久久精品波多野结衣 | 黄色1级毛片 | 亚洲日韩精品欧美一区二区 | 高清在线观看av | 中文字幕在线观看播放 | 欧美一级视频在线观看 | 国产精品永久免费视频 | 在线免费av网站 | 国产黄色片一级 | 成人免费在线看片 | 久久露脸国产精品 | 国产精品嫩草影院9 | 色五月激情五月 | 久久视频在线观看中文字幕 | 久久在草 | 日日日操操 | 久热久草 | 日韩黄色大片在线观看 | 久久日本视频 | 亚洲高清在线观看视频 | 国产精品福利小视频 | 午夜精品电影 | 国产精品四虎 | 又污又黄的网站 | 久久99久久99精品免视看婷婷 | 亚洲激情网站免费观看 | 久影院| 久久综合色影院 | 香蕉网在线观看 | 玖玖爱国产在线 | 中文超碰字幕 | 99精品视频一区二区 | 亚洲天堂色婷婷 | 99精品在这里 | 国产手机在线观看视频 | 国产亚洲精品日韩在线tv黄 | 五月婷婷久久综合 | 顶级bbw搡bbbb搡bbbb | 亚洲综合导航 | 看片一区二区三区 | 国产欧美日韩一区 | 日日干夜夜草 | 日韩在线高清免费视频 | 久久久午夜精品福利内容 | 成人免费视频视频在线观看 免费 | 久久久久成人精品 | 成人黄色电影免费观看 | 日韩三级一区 | 国产精品久久99综合免费观看尤物 | 91精品国产入口 | 成人影视免费 | 午夜精品一区二区三区在线播放 | 亚洲伊人第一页 | 免费试看一区 | 欧美一级片在线观看视频 | 伊人午夜 | 99久久精品无码一区二区毛片 | 在线免费观看黄 | 日日摸日日添日日躁av | 久在线| 欧美动漫一区二区三区 | 久久久久 免费视频 | 97超碰人人澡人人爱 | 国产高清第一页 | 日本久热 | 亚洲影院国产 | 操操操干干干 | 激情欧美日韩一区二区 | 久久撸在线视频 | 国产成人高清av | av大全在线 | 玖玖爱国产在线 | 亚洲片在线资源 | 中文字幕av免费观看 | 青草视频在线 | 涩涩网站免费 | 国产精品网站一区二区三区 | 91免费视频网站在线观看 | 美女视频黄免费网站 | 精品久久久久久久久久 | 欧美日韩首页 | 99色婷婷| 韩国av一区二区 | 91高清完整版在线观看 | 日韩亚洲国产中文字幕 | 久久国产亚洲 | 亚洲伦理精品 | 天天鲁天天干天天射 | 国产精品网红福利 | 碰超在线97人人 | 精品在线一区二区三区 | 4438全国亚洲精品观看视频 | 五月天亚洲婷婷 | 久久96国产精品久久99漫画 | 日韩 在线观看 | 国产精品系列在线 | 小草av在线播放 | 成人黄色小说视频 | 国产精品久久久久久久久久尿 | 狠狠干中文字幕 | 亚洲 中文 欧美 日韩vr 在线 | 国产成人精品午夜在线播放 | 亚洲电影久久久 | 精品久久久久久国产91 | 天天操天天舔天天爽 | 亚洲第一av在线播放 | 色婷婷狠| 久草在线一免费新视频 | 日韩精品国产一区 | 国产精品美女久久久免费 | 日韩黄色一区 | 一区二区在线不卡 | 激情电影影院 | 日韩视频1区| 激情婷婷网 | 高清不卡免费视频 | 中文字幕乱在线伦视频中文字幕乱码在线 | 亚洲欧美日韩一区二区三区在线观看 | 91黄色在线观看 | 国产精品久久三 | 国产精品久久久久9999吃药 | 日日夜夜精品视频 | 韩国av在线播放 | 成人免费在线视频观看 | 黄网站免费大全入口 | 99视频国产精品 | 国产v在线播放 | 欧美成人在线免费观看 | 婷婷 中文字幕 | 一区二区三区在线观看免费视频 | 国产精品久久一区二区三区, | 欧美性另类 | 精品国产一区二区三区av性色 | 日韩免费一二三区 | 亚洲美女精品区人人人人 | 成年人免费av网站 | www黄免费 | www.色午夜| 亚洲精品在线视频 | 不卡的av中文字幕 | 伊人色播| 亚洲视频1区2区 | 六月丁香社区 | 久久大视频 | 狠狠久久综合 | 免费网站黄 | 久久精品久久精品 | 福利二区视频 | 在线观看中文字幕网站 | 久久99精品视频 | 成人国产电影在线观看 | 亚洲精品www久久久久久 | 久久在线观看 | 久草视频在线免费看 | 亚洲精选久久 | 欧美激情va永久在线播放 | 久草99 | 毛片网站免费在线观看 | 国产亚洲一区二区在线观看 | 国产99中文字幕 | 日韩欧美在线观看 | 绯色av一区 | 综合精品久久久 | 国产精品久久久久久久久久久杏吧 | 在线观看免费av网 | 国产无区一区二区三麻豆 | 国产午夜精品一区二区三区在线观看 | 中文字幕 国产精品 | 精品国产观看 | 国模精品一区二区三区 | 成人播放器 | 日本最新一区二区三区 | 日韩欧美一区二区在线观看 | 国产亚洲精品久久久久久大师 | 国产高h视频 | 亚洲午夜激情网 | 久久精品成人 | 国产精品青青 | 久久99九九99精品 | 国产最新视频在线 | 国产色在线 | 久久久久免费 | 97人人爽人人 | 六月久久婷婷 | 亚洲精品高清在线观看 | 日日夜夜网站 | 亚洲精品免费观看 | 天天操天 | 色人久久 | 久久精品视频在线观看 | 91精品国产乱码久久 | 国产精品日韩在线 | 久久久久久久久久久电影 | www.99久久.com | 丁香导航| 伊人久久国产精品 | 欧美色图一区 | 久久99免费| 国产精品高清在线观看 | 亚洲精品456在线播放第一页 | 能在线观看的日韩av | 久久首页| 国产精品wwwwww | 国产精品免费在线播放 | 激情导航| 国产字幕在线观看 | 欧美日韩在线网站 | 中文字幕在线中文 | 精品一区 精品二区 | 欧美日韩国产区 | 久久久午夜视频 | 亚洲经典中文字幕 | 91精品视频免费看 | 国产糖心vlog在线观看 | 国产日韩一区在线 | 色吊丝在线永久观看最新版本 | 精品国产伦一区二区三区观看说明 | 91porny九色91啦中文 | 一区二区三区免费在线观看 | www99久久 | 性色av免费在线观看 | 久久免费毛片 | 91av在线视频播放 | 精品国产亚洲日本 | 久久久精品免费看 | 色欧美88888久久久久久影院 | 97在线成人 | 久久毛片高清国产 | 日韩av播放在线 | 久久久久久久久久影院 | 大胆欧美gogo免费视频一二区 | 国内精品美女在线观看 | 日日干,天天干 | 黄色av观看 | 美女福利视频网 | aⅴ精品av导航 | 国产资源中文字幕 | 色资源中文字幕 | 美女激情影院 | 在线观看成年人 | 成 人 黄 色 视频播放1 | 久草网站 | 91九色网站 | 日韩一级片观看 | 狠狠干综合 | 欧美成人区 | 日韩在线三级 | 超碰人人在 | 又爽又黄又刺激的视频 | 免费高清看电视网站 | 91av手机在线 | 精品国产99 | 国产一级片一区二区三区 | 少妇bbr搡bbb搡bbb | 九九免费在线视频 | 久久国产午夜精品理论片最新版本 | 国产成人在线综合 | 在线播放亚洲 | 在线观看国产中文字幕 | 亚洲成人高清在线 | 日本成人黄色片 | 久久精品视频2 | 91精品国产自产老师啪 | 成人av在线直播 | 天堂网中文在线 | 免费日韩 精品中文字幕视频在线 | 久久国产成人午夜av影院潦草 | 色橹橹欧美在线观看视频高清 | 国产成人精品一区二区三区在线 | 久久精品美女 | 亚洲人人网 | 亚洲欧洲在线视频 | 91在线观看欧美日韩 | 2022久久国产露脸精品国产 | 日韩免费在线视频观看 | 久久成人视屏 | 一区二区中文字幕在线观看 | 91在线看| 久碰视频在线观看 | 免费激情网 | 天天躁天天躁天天躁婷 | 丁香久久激情 | 国产精品99久久久久久久久久久久 | 成人黄色在线观看视频 | 91亚洲国产成人久久精品网站 | 日韩精品视频第一页 | 九九热免费视频在线观看 | 97视频在线观看免费 | 色多多视频在线观看 | 久草网站在线 | 天天操天天操天天干 | 91一区啪爱嗯打偷拍欧美 | 中文在线字幕免 | 日韩av成人在线 | 中文乱幕日产无线码1区 | 精品一二三区 | 色天天中文| 精品久久国产 | 欧美一区二区三区不卡 | 日韩精品一区二区三区免费视频观看 | 午夜精品福利一区二区 | 麻豆系列在线观看 | 国产免费成人av | 精品在线观看视频 | 美女在线观看av | 999超碰| 韩国精品福利一区二区三区 | 国产精品久久久久久久久久久久午夜片 | 成年美女黄网站色大片免费看 | 久久中文欧美 | 国产 日韩 在线 亚洲 字幕 中文 | 中文在线中文a | 黄在线 | 日韩视频三区 | 丁香 婷婷 激情 | 亚洲日韩中文字幕 | 亚洲九九 | 欧美在线视频一区二区三区 | 99精品电影| 久久久一本精品99久久精品 | 免费黄色网址大全 | 亚洲草视频 | 午夜精品久久久久久久99水蜜桃 | 成人av电影免费在线播放 | 91网免费观看 | 国产精品综合av一区二区国产馆 | 最新av在线免费观看 | 亚洲日本欧美在线 | 日色在线视频 | 一区二区不卡视频在线观看 | 久久免费国产电影 | 久久99久久99久久 | 51久久夜色精品国产麻豆 | 黄色高清视频在线观看 | 国产视频综合在线 | 韩国av电影在线观看 | 天天操天天干天天插 | 麻豆视频在线免费观看 | 国产日韩三级 | 日韩免费电影 | 久草视频一区 | 91在线精品秘密一区二区 | 亚洲欧洲视频 | 欧美aaaxxxx做受视频 | 国产在线精品播放 | 亚洲视频网站在线观看 | 久久久久激情视频 | 国产一区二区在线观看视频 | 久草在线在线精品观看 | 国产成年人av | 99九九99九九九视频精品 | 国产在线精品国自产拍影院 | 日韩理论在线 | 日韩在线观看 | 91成人免费电影 | 中文字幕在线观 | 久草影视在线 | 久久视频国产 | 色视频一区 | 欧美日本在线视频 | 97超碰人人澡人人 | 色综合久久88色综合天天 | 日韩在线免费电影 | 国产香蕉久久精品综合网 | 国产美女网站在线观看 | 成年人免费看 | 久久免费精品视频 | 天天激情综合 | 久久精品站 | 天天射天天色天天干 | 成人免费看片网址 | 午夜精品一区二区三区免费 | 香蕉蜜桃视频 | 久久经典视频 | 国产精品免费一区二区三区 | 91成人免费观看视频 | 欧美日本一二三 | 免费观看性生活大片3 | 亚洲a成人v | 超碰com | 久插视频| 中文字幕第一页在线vr | 国产精品热 | 婷婷深爱网 | 国产精品久久久久9999吃药 | www激情网 | 九九免费观看视频 | 国产精品欧美一区二区 | a在线播放| 国产99久久久精品 | 天天鲁一鲁摸一摸爽一爽 | 在线视频 亚洲 | 中文字幕 国产专区 | 欧美日韩在线第一页 | 日本精品在线 | 国产成人久久77777精品 | 在线观看国产麻豆 | 国内一级片在线观看 | 久久免费激情视频 | 91网免费观看 | 一区二区三区高清 | 中文字幕在线久一本久 | 亚洲爱视频 | 国产精品日韩 | 亚洲天天干 | 国产精品久久久久久高潮 | 免费高清影视 | 我爱av激情网 | 亚洲第一区在线观看 | 黄色小说网站在线 | 亚洲国产精品视频在线观看 | 91精品国 | 日韩网站在线免费观看 | 午夜久久久久久久久久影院 | 精品亚洲视频在线 | 久久婷婷一区二区三区 | 久久久久一区二区三区四区 | www.av在线.com | 国产成人久久 | 欧美老女人xx | 欧美另类xxxx | 国产成人精品一区二区三区在线观看 | 国产视频69 |