切换到宽版
阿果石油APP正式发布 石油推广咨询 石油地质QQ群：181983917 物探地震反演群：165218480
阿果石油官方微信公众号：agoil2005      
发帖 回复
返回列表
上一页1...456789下一页
到第
页
  • 9098阅读
  • 83回复

【 1,000 元奖励】If You Firstly 绘出渗流力学教材上横空出世的注采舌进图 [复制链接]

 上一主题 下一主题
离线凤非烟
果园中坚
发帖
1035
财富
1997
威望
6
交易币
0
只看该作者 59楼 发表于: 2013-12-30 | 石油求职招聘就上: 阿果石油英才网
会编程的赶紧上啊，支持一下
回复
举报
离线dragon2009
果园中坚
发帖
2499
财富
436
威望
11
交易币
0
只看该作者 60楼 发表于: 2013-12-30 | 石油求职招聘就上: 阿果石油英才网

加一个线条比较密集的图：

《程序源代码》
clc;
clear;
hold off
%此程序用于绘制舌进现象示意图
%% 基本参数
a=100;      %井距一半,m
q=30;       %产量
%% 绘制坐标轴
%y轴
x=zeros(1,100);
y=linspace(-2*a,2*a,100);
plot(x,y,'k')
hold on
% x轴
y=zeros(1,100);
x=linspace(-3*a,3*a,100);
plot(x,y,'k')
hold on
%% 绘制油水井点
rw=5;   %井点坐标大小
n=50;  %半径劈分次数
m=50;  %角度劈分次数
dr=rw/n;
ds=2*3.14/m;
for i=1:n  
    for j=1:m  
        xw1(i,j)=-a+dr*(i-1)*cos(ds*(j-1));
        yw1(i,j)=dr*(i-1)*sin(ds*(j-1));
        xw2(i,j)=a+dr*(i-1)*cos(ds*(j-1));
        yw2(i,j)=dr*(i-1)*sin(ds*(j-1));
    end
    plot(xw1(i,:),yw1(i,:),'r')
    hold on
    plot(xw2(i,:),yw2(i,:),'b')
    hold on
end
%% 绘制舌进曲线
m2=100;  %角度
ds=2*3.14/m2;
nt=200;
dt=20;
r1(1,1:m2+1)=rw;
flag=0;
for k=2:10000
    for j=1:m2+1
        r2(1,j)=sqrt(r1(1,j)*r1(1,j)+2*a*2*a-2*2*a*r1(1,j)*cos(ds*(j-1)));
        v(k,j)=q*a/(3.14*r1(k-1,j)*r2(k-1,j));
        r1(k,j)=r1(k-1,j)+v(k,j)*dt;
        if r1(k,1)>200
            flag=1;
            T=(k-1)*dt;
            K=k-1;  %次数
        end
       if flag==0
        r2(k,j)=sqrt(r1(k,j)*r1(k,j)+2*a*2*a-2*2*a*r1(k,j)*cos(ds*(j-1)));
       else
           break
       end
    end
    if flag==1
        break
    end
end
%计算个点坐标
for k=1:K
    for j=1:m2+1
        x(k,j)=-a+r1(k,j)*cos(ds*(j-1));
        y(k,j)=r1(k,j)*sin(ds*(j-1));
    end
    plot(x(k,:),y(k,:),'m')
    hold on
end
学石油、爱石油，我为祖国献石油……
回复
举报
离线齊成偉
专栏作者
发帖
663
财富
42
威望
28
交易币
0
只看该作者 61楼 发表于: 2013-12-30 | 石油求职招聘就上: 阿果石油英才网
回 dragon2009 的帖子
dragon2009:[图片]
加一个线条比较密集的图：
[图片]
 (2013-12-30 21:41) 

加我QQ 1261621
回复
举报
离线油猫警长
果园新丁
发帖
52
财富
30
威望
0
交易币
0
只看该作者 62楼 发表于: 2013-12-30 | 石油求职招聘就上: 阿果石油英才网
代码：


%一部分流线
clc;clear
d = 10;
hold on
for i = 1:100
c0 = 0.01*i;
c = (1-c0)/(2*c0);
theta0 = atan(1/c);
rs = d^2*(1+1/c^2);
r = sqrt(rs);
theta = -theta0:0.1:pi+theta0;
%theta = 0:0.01:2*pi;
x = r*cos(theta);
y = d/c+r*sin(theta);
plot(x,y)
plot(x,-y)
end
grid on
axis equal
axis([-35 35 -35 35])
hold off
%另一部分流线
clc;clear
d = 10;
hold on
for i = 1:100
c0 = 0.01*i;
c = (1-c0)/(2*c0);
theta0 = atan(1/c);
rs = d^2*(1+1/c^2);
r = sqrt(rs);
theta = theta0:0.1:pi-theta0;
%theta = 0:0.01:2*pi;
x = r*cos(theta);
y = -d/c+r*sin(theta);
plot(x,y)
plot(x,-y)
end
grid on
axis equal
axis([-35 35 -35 35])
hold off
%绘制舌进
clc
clear
al = 1;
d = 10;
hold on
t0 = 500;
for i = 1:70
c0 = 0.01*i;
c = (1-c0)/(2*c0);
theta0 = atan(1/c);
rs = d^2*(1+1/c^2);
r = sqrt(rs);
theta1 = -theta0:0.01:pi;
n = length(theta1);
for i = 1:n
    t = quad(@(theta)fun(theta,r,d,c,al),-theta0,theta1(i));
    if t >= t0
        the = theta1(i);
        break
    end
end
x = r*cos(the);
y = r*sin(the)+d/c;
plot(x,y,'r.')
plot(x,-y,'r.')
end

for i = 71:98
c0 = 0.01*i;
c = (1-c0)/(2*c0);
theta0 = atan(1/c);
rs = d^2*(1+1/c^2);
r = sqrt(rs);
theta1 = -theta0:0.0001:pi;
n = length(theta1);
for i = 1:n
    t = quad(@(theta)fun(theta,r,d,c,al),-theta0,theta1(i));
    if t >= t0
        the = theta1(i);
        break
    end
end
x = r*cos(the);
y = r*sin(the)+d/c;
plot(x,y,'r.')
plot(x,-y,'r.')
end
%2
for i = 1:100
c0 = 0.01*i;
c = (1-c0)/(2*c0);
theta0 = atan(1/c);
rs = d^2*(1+1/c^2);
r = sqrt(rs);
theta2 = theta0:0.01:pi;
n = length(theta2);
for i = 1:n
    t = quad(@(theta)fun2(theta,r,d,c,al),theta0,theta2(i));
    if t >= t0
        the = theta2(i);
        break
    end
end
x = r*cos(the);
y = r*sin(the)-d/c;
plot(x,y,'r.')
plot(x,-y,'r.')
end

%调用函数1
function f = fun(theta,r,d,c,al)
f = r/al*(sqrt((r*cos(theta)-d).^2+(r*sin(theta)+d/c).^2).*sqrt((r*cos
(theta)+d).^2+(r*sin(theta)+d/c).^2));
%调用函数2
function f = fun2(theta,r,d,c,al)
f = r/al*(sqrt((r*cos(theta)-d).^2+(r*sin(theta)-d/c).^2).*sqrt((r*cos
(theta)+d).^2+(r*sin(theta)-d/c).^2));
回复
举报
离线齊成偉
专栏作者
发帖
663
财富
42
威望
28
交易币
0
只看该作者 63楼 发表于: 2013-12-30 | 石油求职招聘就上: 阿果石油英才网
回 油猫警长 的帖子
油猫警长:代码：
%一部分流线
clc;clear
....... (2013-12-30 22:49)

【复制备份 以防更改】我这个裁判不好当 请谅解
代码：


%一部分流线
clc;clear
d = 10;
hold on
for i = 1:100
c0 = 0.01*i;
c = (1-c0)/(2*c0);
theta0 = atan(1/c);
rs = d^2*(1+1/c^2);
r = sqrt(rs);
theta = -theta0:0.1:pi+theta0;
%theta = 0:0.01:2*pi;
x = r*cos(theta);
y = d/c+r*sin(theta);
plot(x,y)
plot(x,-y)
end
grid on
axis equal
axis([-35 35 -35 35])
hold off
%另一部分流线
clc;clear
d = 10;
hold on
for i = 1:100
c0 = 0.01*i;
c = (1-c0)/(2*c0);
theta0 = atan(1/c);
rs = d^2*(1+1/c^2);
r = sqrt(rs);
theta = theta0:0.1:pi-theta0;
%theta = 0:0.01:2*pi;
x = r*cos(theta);
y = -d/c+r*sin(theta);
plot(x,y)
plot(x,-y)
end
grid on
axis equal
axis([-35 35 -35 35])
hold off
%绘制舌进
clc
clear
al = 1;
d = 10;
hold on
t0 = 500;
for i = 1:70
c0 = 0.01*i;
c = (1-c0)/(2*c0);
theta0 = atan(1/c);
rs = d^2*(1+1/c^2);
r = sqrt(rs);
theta1 = -theta0:0.01:pi;
n = length(theta1);
for i = 1:n
    t = quad(@(theta)fun(theta,r,d,c,al),-theta0,theta1(i));
    if t >= t0
        the = theta1(i);
        break
    end
end
x = r*cos(the);
y = r*sin(the)+d/c;
plot(x,y,'r.')
plot(x,-y,'r.')
end

for i = 71:98
c0 = 0.01*i;
c = (1-c0)/(2*c0);
theta0 = atan(1/c);
rs = d^2*(1+1/c^2);
r = sqrt(rs);
theta1 = -theta0:0.0001:pi;
n = length(theta1);
for i = 1:n
    t = quad(@(theta)fun(theta,r,d,c,al),-theta0,theta1(i));
    if t >= t0
        the = theta1(i);
        break
    end
end
x = r*cos(the);
y = r*sin(the)+d/c;
plot(x,y,'r.')
plot(x,-y,'r.')
end
%2
for i = 1:100
c0 = 0.01*i;
c = (1-c0)/(2*c0);
theta0 = atan(1/c);
rs = d^2*(1+1/c^2);
r = sqrt(rs);
theta2 = theta0:0.01:pi;
n = length(theta2);
for i = 1:n
    t = quad(@(theta)fun2(theta,r,d,c,al),theta0,theta2(i));
    if t >= t0
        the = theta2(i);
        break
    end
end
x = r*cos(the);
y = r*sin(the)-d/c;
plot(x,y,'r.')
plot(x,-y,'r.')
end

%调用函数1
function f = fun(theta,r,d,c,al)
f = r/al*(sqrt((r*cos(theta)-d).^2+(r*sin(theta)+d/c).^2).*sqrt((r*cos
(theta)+d).^2+(r*sin(theta)+d/c).^2));
%调用函数2
function f = fun2(theta,r,d,c,al)
f = r/al*(sqrt((r*cos(theta)-d).^2+(r*sin(theta)-d/c).^2).*sqrt((r*cos
(theta)+d).^2+(r*sin(theta)-d/c).^2));

【复制备份 以防更改】我这个裁判不好当 请谅解
回复
举报
离线齊成偉
专栏作者
发帖
663
财富
42
威望
28
交易币
0
只看该作者 64楼 发表于: 2013-12-30 | 石油求职招聘就上: 阿果石油英才网
回 dragon2009 的帖子
dragon2009:
[图片]
加一个线条比较密集的图：
[图片]
《程序源代码》
clc;
.......

【复制备份 以防更改】我这个裁判不好当 请谅解

《程序源代码》clc;
clear;
hold off
%此程序用于绘制舌进现象示意图
%% 基本参数
a=100;      %井距一半,m
q=30;       %产量
%% 绘制坐标轴
%y轴
x=zeros(1,100);
y=linspace(-2*a,2*a,100);
plot(x,y,'k')
hold on
% x轴
y=zeros(1,100);
x=linspace(-3*a,3*a,100);
plot(x,y,'k')
hold on
%% 绘制油水井点
rw=5;   %井点坐标大小
n=50;  %半径劈分次数
m=50;  %角度劈分次数
dr=rw/n;
ds=2*3.14/m;
for i=1:n  
    for j=1:m  
        xw1(i,j)=-a+dr*(i-1)*cos(ds*(j-1));
        yw1(i,j)=dr*(i-1)*sin(ds*(j-1));
        xw2(i,j)=a+dr*(i-1)*cos(ds*(j-1));
        yw2(i,j)=dr*(i-1)*sin(ds*(j-1));
    end
    plot(xw1(i,:),yw1(i,:),'r')
    hold on
    plot(xw2(i,:),yw2(i,:),'b')
    hold on
end
%% 绘制舌进曲线
m2=100;  %角度
ds=2*3.14/m2;
nt=200;
dt=20;
r1(1,1:m2+1)=rw;
flag=0;
for k=2:10000
    for j=1:m2+1
        r2(1,j)=sqrt(r1(1,j)*r1(1,j)+2*a*2*a-2*2*a*r1(1,j)*cos(ds*(j-1)));
        v(k,j)=q*a/(3.14*r1(k-1,j)*r2(k-1,j));
        r1(k,j)=r1(k-1,j)+v(k,j)*dt;
        if r1(k,1)>200
            flag=1;
            T=(k-1)*dt;
            K=k-1;  %次数
        end
       if flag==0
        r2(k,j)=sqrt(r1(k,j)*r1(k,j)+2*a*2*a-2*2*a*r1(k,j)*cos(ds*(j-1)));
       else
           break
       end
    end
    if flag==1
        break
    end
end
%计算个点坐标
for k=1:K
    for j=1:m2+1
        x(k,j)=-a+r1(k,j)*cos(ds*(j-1));
        y(k,j)=r1(k,j)*sin(ds*(j-1));
    end
    plot(x(k,:),y(k,:),'m')
    hold on
end

【复制备份 以防更改】我这个裁判不好当 请谅解
回复
举报
离线油猫警长
果园新丁
发帖
52
财富
30
威望
0
交易币
0
只看该作者 65楼 发表于: 2013-12-30 | 石油求职招聘就上: 阿果石油英才网


以上是公式推导
回复
举报
离线齊成偉
专栏作者
发帖
663
财富
42
威望
28
交易币
0
只看该作者 66楼 发表于: 2013-12-31 | 石油求职招聘就上: 阿果石油英才网
回 油猫警长 的帖子

油猫警长:
[图片]
以上是公式推导

【复制备份 以防更改】我这个裁判不好当 请谅解








【复制备份 以防更改】我这个裁判不好当 请谅解
回复
举报
离线齊成偉
专栏作者
发帖
663
财富
42
威望
28
交易币
0
只看该作者 67楼 发表于: 2013-12-31 | 石油求职招聘就上: 阿果石油英才网
回 dragon2009 的帖子
dragon2009:
[图片]
加一个线条比较密集的图：
[图片]
《程序源代码》
clc;
.......

【复制备份 以防更改】我这个裁判不好当 请谅解


【复制备份 以防更改】我这个裁判不好当 请谅解
回复
举报
离线newlang
果园中坚
发帖
1057
财富
2525
威望
21
交易币
0
只看该作者 68楼 发表于: 2013-12-31 | 石油求职招聘就上: 阿果石油英才网
谢谢鼓励
回复
举报
发帖 回复
返回列表
上一页1...456789下一页
到第
页

网站事务咨询:QQ：1392013 | 26189883
阿果石油网为免费个人网站，为石油人提供免费的在线即时技术交流场所，拒绝任何人以任何形式在本论坛发表与中华人民共和国法律相抵触的言论和行为!
如有言论或会员共享的资料涉及到您的权益，请立即通知网站管理员，本站将在第一时间给予配合处理，谢谢！