欢迎访问 生活随笔!

生活随笔

当前位置: 首页 > 编程资源 > 编程问答 >内容正文

编程问答

潮汐特征值计算及其预报

发布时间:2023/12/20 编程问答 38 豆豆
生活随笔 收集整理的这篇文章主要介绍了 潮汐特征值计算及其预报 小编觉得挺不错的,现在分享给大家,帮大家做个参考.

        本人为海洋大学大二在读学生,由于作业需求。自行编写matlab程序以计算潮汐类型书并对潮汐进行了预报。

        我国作为陆地大国,由于各种历史问题,在过去忽视了海洋及海洋产业的发展,这也导致目前我国海洋学正处于起步阶段,网上各类海洋学相关资料也是少之又少。对于海洋学子,没有外部资料无疑增加了海洋学的学习难度。而由于本人并非计算机相关专业,并未受过相关数据结构,编程风格等计算机的专业训练,因此本人代码可能会有不规范的地方。

        并且,本人计算出的结果有较大误差,由于本人能力有限,无法找到错误所在。希望本篇文章抛砖引玉,能够引起大家的思考和讨论,不断完善代码,改进算法。希望大家多加注释,未后继的海洋学生提供学习的素材。

以下为程序的主函数部分:(潮汐观测文件会在后续上传)

close all; clear;%% 导入潮汐数据 % 第一个数对应的时间为2003年3月3日0时,相邻数据的时间间隔为1个小时。 tide_arr = importdata('(本地潮汐数据)'); %此处为导入本地潮汐数据文件%% 根据已知数据绘出潮流观测曲线 x1 = 1:721; %x为观测矩阵中第一列,既数据初始编号 y1 = tide_arr(:,2); %y为观测矩阵中第二列,既观测潮位 plot(x1,y1); xlabel('序列号'); ylabel('潮位'); title('观测潮位与自报潮位'); hold on %% 选取8个主要分潮O1,P1,K1,M2,K2,Q1,S2,N2 % 查表得八个主要分潮的杜德森数 % 数组中分别对应u1,u2,u3,u4,u5,u6,u0 O1_du = [1,-1,0,0,0,0,-1]; P1_du = [1,1,-2,0,0,0,-1]; K1_du = [1,1,0,0,0,0,1]; M2_du = [2,0,0,0,0,0,0]; K2_du = [2,2,0,0,0,0,0]; Q1_du = [1,-2,0,1,0,0,-1]; S2_du = [2,2,-2,0,0,0,0]; N2_du = [2,-1,0,1,0,0,0]; %% 计算时间原点 % 第一个数对应的时间为2003年3月3日0时,相邻数据的时间间隔为1个小时。 year=2003; month=3; day=3; hour=0; time_interval = 1; time_calculator(year,month,day,hour,time_interval); %计算出中间时刻对应的时间并输出 %% 计算基本天文元素随时间的变化速度 omega_t =(2*pi) / (24 + (50/60)); %平太阴时角,周期一个平太阴日=24小时50分钟 omega_s = (2*pi) / (27.32*24); %月球平均经度,周期一个回归月 = 27.32个平太阳日 omega_h = (2*pi) / (365.2422*24); %太阳平均经度,周期一个回归年 = 265.2422个平太阳日 omega_p1 = (2*pi) / (8.85*365.2422*24); %月球近地点平均经度,周期为8.85年 omega_n = (2*pi) / (18.61*365.2422*24); %白道升交点在黄道上经度的负值,周期为18.61年 omega_p2 = (2*pi) / (21000*365.2422*24); %太阳近地点平均经度,周期为21000年,u6=0此项可以忽略%基本天文分潮随时间变化率矩阵 basic_value_omega = [omega_t,omega_s,omega_h,omega_p1,omega_n,omega_p2]'; %% 分别计算八个分潮的角速率 omega_O1 = O1_du(1:6) * basic_value_omega; omega_P1 = P1_du(1:6) * basic_value_omega; omega_K1 = K1_du(1:6) * basic_value_omega; omega_M2 = M2_du(1:6) * basic_value_omega; omega_K2 = K2_du(1:6) * basic_value_omega; omega_Q1 = Q1_du(1:6) * basic_value_omega; omega_S2 = S2_du(1:6) * basic_value_omega; omega_N2 = N2_du(1:6) * basic_value_omega; %创建八个分潮角速度向量 %Omega = [omega_O1,omega_P1,omega_K1,omega_M2,omega_K2,omega_Q1,omega_S2,omega_N2];% fprintf('O1分潮的角速率:%frad/h\n',omega_O1); % fprintf('P1分潮的角速率:%frad/h\n',omega_P1); % fprintf('K1分潮的角速率:%frad/h\n',omega_K1); % fprintf('M2分潮的角速率:%frad/h\n',omega_M2); % fprintf('K2分潮的角速率:%frad/h\n',omega_K2); % fprintf('Q1分潮的角速率:%frad/h\n',omega_Q1); % fprintf('S2分潮的角速率:%frad/h\n',omega_S2); % fprintf('N2分潮的角速率:%frad/h\n',omega_N2);% 从这往上计算均正确

        代码从这以上应该没有问题,以下的初相计算和标准答案有较大偏差,本人能力有限,对知识的理解不够,未能找出错误。

%% 计算天文元素(以弧度制计算) %初项为第一个数据收集时刻的位相:2003年3月3日0时 Y = 2003; M = 3; D = 3; t = 0; n = D + 58; i = fix((Y - 1900)/4); %月球平均经度s s = (pi/180)*(277.02) + (pi/180)*(129.3848)*(Y-1900) + (pi/180)*(13.1764)*(n+i+t/24); %太阳平均经度h h = (pi/180)*(280.19) + (pi/180)*(0.2387)*(Y-1900) + (pi/180)*(0.9857)*(n+i+t/24); %月球近地点平均经度p1 p1 = (pi/180)*(334.39) + (pi/180)*(40.6625)*(Y-1900) + (pi/180)*(0.1114)*(n+i+t/24); %白道升交点在黄道上经度的负值n n = (pi/180)*(100.84) + (pi/180)*(19.3282)*(Y-1900) + (pi/180)*(0.0530)*(n+i+t/24); %太阳近地点平均经度p2 p2 = (pi/180)*(281.22) + (pi/180)*(0.0172)*(Y-1900) + (pi/180)*(0.00005)*(n+i+t/24); %τ=15t?s+h′ 平太阴时tau tau=(pi/180)*15*t-s+h; fprintf('平太阴时角τ:%f\n月球平均经度s:%f\n太阳平均经度h:%f\n月球近地点平均经度p1:%f\n白道升交点在黄道上经度的负值n:%f\n太阳近地点平均经度p2:%f\n',tau,s,h,p1,n,p2); basic_value = [tau,s,h,p1,n,p2,pi/2]'; %基本天文元素矩阵% 调整到[0,2*pi] %将初始位相调整到[0,2*pi] for i = 1:7if basic_value(i) > 2*piwhile basic_value(i) > 2*pibasic_value(i) = basic_value(i)-2*pi;endelseif basic_value(i) < (-2)*piwhile basic_value(i) < (-2)*pibasic_value(i) = basic_value(i)+2*pi;endendif basic_value(i) < 0 basic_value(i) = basic_value(i) + 2*pi;end end% s=14934.4708; s=deg2rad(s); % h=355.1596; h=deg2rad(h); % p1=4533.8789; p1=deg2rad(p1); % n=2096.9976; n=deg2rad(n); % p2=282.99665; p2=deg2rad(p2); % tau=(pi/180)*15*0-s+h; % basic_value = [tau,s,h,p1,n,p2,pi/2]'; %% 计算初始位相 O1_v0 = O1_du * basic_value; P1_v0 = P1_du * basic_value; K1_v0 = K1_du * basic_value; M2_v0 = M2_du * basic_value; K2_v0 = K2_du * basic_value; Q1_v0 = Q1_du * basic_value; S2_v0 = S2_du * basic_value; N2_v0 = N2_du * basic_value; % fprintf('O1分潮初始位相:%f\n',O1_v0); % fprintf('P1分潮初始位相:%f\n',P1_v0); % fprintf('K1分潮初始位相:%f\n',K1_v0); % fprintf('M2分潮初始位相:%f\n',M2_v0); % fprintf('K2分潮初始位相:%f\n',K2_v0); % fprintf('Q1分潮初始位相:%f\n',Q1_v0); % fprintf('S2分潮初始位相:%f\n',S2_v0); % fprintf('N2分潮初始位相:%f\n',N2_v0);%建立初项的矩阵 all_v0 = [O1_v0,P1_v0,K1_v0,M2_v0,K2_v0,Q1_v0,S2_v0,N2_v0];%将初始位相调整到[0,2*pi] for i = 1:8if all_v0(i) > 2*piwhile all_v0(i) > 2*piall_v0(i) = all_v0(i)-2*pi;endelseif all_v0(i) < (-2)*piwhile all_v0(i) < (-2)*piall_v0(i) = all_v0(i)+2*pi;endendif all_v0(i) < 0 all_v0(i) = all_v0(i) + 2*pi;end end%all_v0 = [0.0240 0.0000 0.7119 6.1142 1.4863 4.8209 4.7969 5.5087]; %% 计算焦点订正角 %基本方法计算方法 %K1分潮的焦点因子f_K1和焦点订正角u_K1 a = [p1,n,p2,pi/2]'; K1_fcosu = 0.0002*cos([-2,-1,0,0]*a) + 0.0001*cos([0,-2,0,0]*a)....+0.0198*cos([0,-1,0,-2]*a) + 1*cos([0,0,0,0]*a)...+0.1356*cos([0,1,0,0]*a) + 0.0029*cos([0,2,0,-2]*a);K1_fsinu = 0.0002*sin([-2,-1,0,0]*a) + 0.0001*sin([0,-2,0,0]*a)....+0.0198*sin([0,-1,0,-2]*a) + 1*sin([0,0,0,0]*a)...+0.1356*sin([0,1,0,0]*a) + 0.0029*sin([0,2,0,-2]*a);f_K1 = sqrt( (K1_fcosu)^2 + (K1_fsinu)^2 ); u_K1 = atan(K1_fsinu/K1_fcosu); fprintf('K1分潮的焦点因子:%f\n和焦点订正角:%f\n',f_K1,u_K1); %除K1分潮其余分潮通过自定义函数f_and_u_calculator计算%% 用自定义函数计算其余分潮的焦点因子f和焦点订正角u %O1分潮 p_O1 = [-0.0058,0.1885,1,0.0002,-0.0064,-0.0010]; %O1分潮引潮力系数比 delta_u_01 = [0,-2,0,0;0,-1,0,0;0,0,0,0;2,-1,0,0;2,0,0,0;2,1,0,0]; %O1分潮次要分潮与主要分潮的差值Δμ fprintf('O1分潮的'); [f_O1,u_O1] = f_and_u_calculator(p_O1,delta_u_01);%P1分潮 p_P1 = [0.0008,-0.0112,1,-0.0015,-0.0003]; %P1分潮引潮力系数比 delta_u_P1 = [0,-2,0,0;0,-1,0,0;0,0,0,0;2,0,0,0;2,1,0,0]; %P1分潮次要分潮与主要分潮的差值Δμ fprintf('P1分潮的'); [f_P1,u_P1]= f_and_u_calculator(p_P1,delta_u_P1);%M2分潮 p_M2 = [0.0005,-0.0373,1,0.0006,0.0002]; %M2分潮引潮力系数比 delta_u_M2 = [0,-2,0,0;0,-1,0,0;0,0,0,0;2,0,0,0;2,1,0,0]; %M2分潮次要分潮与主要分潮的差值Δμ fprintf('M2分潮的'); [f_M2,u_M2] = f_and_u_calculator(p_M2,delta_u_M2);%K2分潮 p_K2 = [-0.0128,1,0.2980,0.0324]; %K2分潮引潮力系数比 delta_u_K2 = [0,-1,0,0;0,0,0,0;0,1,0,0;0,2,0,0]; %K2分潮次要分潮与主要分潮的差值Δμ fprintf('K2分潮的'); [f_K2,u_K2] = f_and_u_calculator(p_K2,delta_u_K2);%% Q1,S2,N2三个分潮通过查表得到 %Q1分潮 %Q1分潮的焦点因子和焦点订正角都和O1分潮相同 fprintf('Q1分潮的'); [f_Q1,u_Q1] = f_and_u_calculator(p_O1,delta_u_01); %与O1分潮的相同%S2分潮 %S2分潮的焦点因子为1,焦点订正角为0 f_S2 = 1; u_S2 = 0; fprintf('S2分潮的焦点因子为:1\n焦点订正角为:0\n');%N2分潮 %N2分潮的焦点因子和焦点订正角都和M2分潮相同 fprintf('N2分潮的'); [f_N2,u_N2] = f_and_u_calculator(p_M2,delta_u_M2); %与O1分潮的相同% 焦点因子矩阵 all_f = [f_O1,f_P1,f_K1,f_M2,f_K2,f_Q1,f_S2,f_N2]; %焦点订正角矩阵 all_u = [u_O1,u_P1,u_K1,u_M2,u_K2,u_Q1,u_S2,u_N2]; %% 联立矛盾方程,利用最小二乘法得到法方程并求出矛盾方程最优解 %一共选取8个分潮 %一共有721个观测记录数据,那么就有721个方程构成矛盾方程,观测水位h有721个 %每个分潮的振幅和迟角为调和常数,一单确定则不随时间变化,平均水位也未知 %因此未知量的个数有2×8+1=17个(8个分潮的Hsinθ和Hcosθ以及S0) %直接求法方程的系数矩阵A,B,F′,F"%创建八个分潮角速度向量 N = 721 ; Omega = [omega_O1,omega_P1,omega_K1,omega_M2,omega_K2,omega_Q1,omega_S2,omega_N2]; %此步需要注意后面分潮计算的顺序都要按照上面的顺序%% 下面的循环计算系数矩阵A A = zeros(8,8); %初始化 %A矩阵的第一行第一个 A(1,1) = N; %A矩阵除去第一行和第一列 for i = 2:9for j = 2:9A(j,j) = (1/2) * (N + (sin(N*Omega(j-1)*time_interval)) / (sin(Omega(j-1)*time_interval)) );if i ~= jA(i,j) = (1/2) * ((sin((N/2)*(Omega(i-1)-Omega(j-1))*time_interval)) / ...(sin((1/2)*(Omega(i-1)-Omega(j-1))*time_interval))...+ (sin((N/2)*(Omega(i-1)+Omega(j-1))*time_interval)) /...(sin((1/2)*(Omega(i-1)+Omega(j-1))*time_interval)));end % A(j,i) = A(i,j);end end %计算A矩阵的第一行第二个元素到最后一个元素,以及第一列第二个元素到最后一个元素 for j = 2:9A(1,j) = (sin((N/2)*Omega(j-1)*time_interval)) / (sin((1/2)*Omega(j-1)*time_interval));A(j,1) = A(1,j); end %disp(A) %% 下面的循环计算系数矩阵B B = zeros(8,8); for i = 1:8for j = 1:8B(j,j) = 1/2*(N - (sin(N*Omega(j)*time_interval))/(sin(Omega(j)*time_interval)));if i ~= jB(i,j) = (1/2) * ((sin((N/2)*(Omega(i)-Omega(j))*time_interval)) / ...(sin((1/2)*(Omega(i)-Omega(j))*time_interval))...- (sin((N/2)*(Omega(i)+Omega(j))*time_interval)) /...(sin((1/2)*(Omega(i)+Omega(j))*time_interval)));end % B(j,i) = B(i,j);end end %disp(B)%% 下面计算F'和F",F'在下面用向量F1表示,F"用向量F2表示 tide_arr(:,3)=[-360:0,1:360]; F1 = zeros(1,9)'; F2 = zeros(1,8)'; F1(1) = sum(tide_arr(:,2)); for i = 2:9for n = tide_arr(1,3):tide_arr(end,3)F1(i) = F1(i) + (tide_arr(n+361,2)) * cos(n*Omega(i-1)*time_interval);F2(i-1) = F2(i-1) + (tide_arr(n+361,2))*sin(n*Omega(i-1)*time_interval);end end%计算未知量X=Hcosθ,Y=Hsinθ X = F1'/A; Y = F2'/B; S0 = X(1); %准调和分潮振幅为 R = sqrt( power(X(:,2:end),2) + power(Y,2) );%t=0时刻的位相 theta = asin(Y ./ R); for i = 1:8if theta(:,i) < 0;theta(:,i) = pi - theta(:,i);end end% 振幅和位相对应的八个分潮分别为 % Omega = [omega_O1,omega_P1,omega_K1,omega_M2,omega_K2,omega_Q1,omega_S2,omega_N2]; % O1,P1,K1,M2,K2,Q1,S2,N2 %% 计算调和常数 % 振幅 H = R ./ all_f; % 迟角 g = theta + all_v0 + all_u;%% 潮汐自报 h = zeros(721,2); for n = tide_arr(:,3)for i = 1:8 h(n+361) = h(n+361) + all_f(i)*H(i)*cos(n*Omega(i)*time_interval + all_v0(i) + all_u(i) - g(i) );end end h = h(:,1) + S0; %加上平均潮位 h(:,2) = 1:721; %给自报潮位编号%自报潮位作图 x2 = h(:,2); y2 = h(:,1); plot(x2,y2); legend('观测','自报');%% 自报余差计算 %余差r r = tide_arr(:,2) - h(:,1); %自报余差的标准差delta_r_2 = (1/721) * sum(power(r,2)); %方差delta_r = sqrt(delta_r_2); %标准差text(0,190,'自报余差的方差为:185.5727 标准差为:13.6225','Color','green','FontSize',14);%% 预报潮位(5月1日0时到6月1日0时) %已有数据第一个数对应的时间为2003年3月3日0时,相邻数据的时间间隔为1个小时。 %721个数据,从第一个观测记录数到最后一个观测记录数一共经历了720个小时 %所以一共经历了30个平太阳日 %已有的观测记录数从2003年3月3日0时到2003年4月3日0时h_forecast = zeros(721,2);% 潮汐的振幅H和迟角g为调和常量,是本身的动力学要素,不随时间改变 % S0一但确定也不会改变 % 焦点因子f和焦点订正角u在1年之内变化缓慢,近似看作不变%计算新的天文元素 %计算待预报潮位的中间时刻,一共也会有721个数据 Y_forecast_middle = 2003 ; M_forecast_middle = 5; t_total = 360; D_forecast_middle = fix(t_total/24); t_forecast_middle = rem(t_total,24); basoc_valie_forecast = zeros(7,1); basic_value_forecast = basic_value_calculator(Y_forecast_middle,M_forecast_middle,D_forecast_middle,t_forecast_middle);%% 计算新的初相位 O1_v0_forecast = O1_du * basic_value_forecast; P1_v0_forecast = P1_du * basic_value_forecast; K1_v0_forecast = K1_du * basic_value_forecast; M2_v0_forecast = M2_du * basic_value_forecast; K2_v0_forecast = K2_du * basic_value_forecast; Q1_v0_forecast = Q1_du * basic_value_forecast; S2_v0_forecast = S2_du * basic_value_forecast; N2_v0_forecast = N2_du * basic_value_forecast; %预报的初相矩阵 all_v0_forecast = [O1_v0_forecast,P1_v0_forecast,K1_v0_forecast,M2_v0_forecast,K2_v0_forecast,Q1_v0_forecast,S2_v0_forecast,N2_v0_forecast];for n = tide_arr(:,3)for i = 1:8 h_forecast(n+361) = h_forecast(n+361) + ...all_f(i)*H(i)*cos(n*Omega(i)*time_interval + all_v0_forecast(i) + all_u(i) - g(i) );end end h_forecast = h_forecast(:,1) + S0; %加上平均潮位 h_forecast(:,2)=1:721;%% 画出预报的潮位 figure(2); x3 = h_forecast(:,2); y3 = h_forecast(:,1); plot(x3,y3); xlabel('预报记录数'); ylabel('潮位'); title('预报潮位(5月1日0时到6月1日0时)');

以下部分为自编函数部分,主要分为三个:

%这个函数用来计算天文元素(以弧度制计算)并返回基本天文元素矩阵 %这个函数未考虑平年和闰年的情况 %记1月1日为第0天function basic_value = basic_value_calculator(Y,M,D,t)%中间时刻的世纪时为:2003年3月18日0时 % Y = 2003; M = 3; D = 18; t = 0; if M == 1n = D - 1; elseif M == 2n = D +30; elseif M == 3n = D + 58; elseif M == 4n = D + 89; elseif M == 5n = D + 119; elseif M == 6n = D + 150; elseif M == 7n = D + 180; elseif M == 8n = D + 211; elseif M == 9n = D + 242; elseif M == 10n = D + 272; elseif M == 11n = D + 303; elseif M == 12n = D + 333; endi = fix((Y - 1900)/4); %月球平均经度s s = (pi/180)*(277.02) + (pi/180)*(129.3848)*(Y-1900) + (pi/180)*(13.1764)*(n+i+t/24); %太阳平均经度h h = (pi/180)*(280.19) + (pi/180)*(0.2387)*(Y-1900) + (pi/180)*(0.9857)*(n+i+t/24); %月球近地点平均经度p1 p1 = (pi/180)*(334.39) + (pi/180)*(40.6625)*(Y-1900) + (pi/180)*(0.1114)*(n+i+t/24); %白道升交点在黄道上经度的负值n n = (pi/180)*(100.84) + (pi/180)*(19.3282)*(Y-1900) + (pi/180)*(0.0530)*(n+i+t/24); %太阳近地点平均经度p2 p2 = (pi/180)*(281.22) + (pi/180)*(0.0172)*(Y-1900) + (pi/180)*(0.00005)*(n+i+t/24); %τ=15t?s+h′ 平太阴时tau tau=(pi/180)*15*t-s+h; %fprintf('平太阴时角τ:%f\n月球平均经度s:%f\n太阳平均经度h:%f\n月球近地点平均经度p1:%f\n白道升交点在黄道上经度的负值n:%f\n太阳近地点平均经度p2:%f\n',tau,s,h,p1,n,p2); basic_value = [tau,s,h,p1,n,p2,pi/2]'; %返回基本天文元素矩阵 %此函数用来计算分潮的焦点因子f和焦点订正角u %函数使用时注意:引潮力系数比和杜德森数之差要一一对应 %ρ前带正负号时,既不考虑μ0时,要将μ0位置用0代替。function [a,b] = f_and_u_calculator(p,delta_u) %设置参数:引潮力系数比为p, a = [p1,n,p2,pi/2]'为一个定值 %delta_u为一个矩阵,是同一亚群其他分潮与主要分潮u4,u5,u6,u0之差,其中u6恒为0 %写delta_u矩阵时要注意每行与引潮力系数比p的分潮要一一对应%% 由于基本天文元素不变,因此将矩阵a设置为定值 Y = 2003; D = 3; t = 0; n = D + 58; i = fix((Y - 1900)/4); p1 = (pi/180)*(334.39) + (pi/180)*(40.6625)*(Y-1900) + (pi/180)*(0.1114)*(n+i+t/24); %白道升交点在黄道上经度的负值n n = (pi/180)*(100.84) + (pi/180)*(19.3282)*(Y-1900) + (pi/180)*(0.0530)*(n+i+t/24); %太阳近地点平均经度p2 p2 = (pi/180)*(281.22) + (pi/180)*(0.0172)*(Y-1900) + (pi/180)*(0.00005)*(n+i+t/24); a = [p1,n,p2,pi/2]';%% 焦点因子和焦点订正角计算的主要步骤 fcosu = 0; fsinu = 0; for j = 1:length(p)fcosu = fcosu + p(j) * cos(delta_u(j,:)*a);fsinu = fsinu + p(j) * sin(delta_u(j,:)*a); endf = sqrt( (fcosu^2)+(fsinu^2) ); u = atan( (fsinu)/(fcosu) );% u=acot(fcosu/fsinu); % f=fcosu/(cos(u)); % u=u+2*pi; % fprintf('焦点因子为:%f\n焦点订正角为:%f\n',f,u); a=f; b=u; end %这个函数用来计算分潮的中间时刻以及其中间时刻的真实时间 %中间时刻N' = 第一个观测记录数 + 1/2 *(N-1)%目前任务的第一个数对应的时间为2003年3月3日0时,相邻数据的时间间隔为1个小时。 %总共有30天的数据function time_calculator(year,month,day,hour,time_interval) tide_arr = importdata('D:\SHOU\Matlab\bin\tide_work\21.txt'); [N,~] = size(tide_arr); %N为总的观测记录数 N_pie = 1+(1/2) * (N-1); %N_pie为中间时刻的编号,第一个观测记录数为1 fprintf('中间时刻潮汐编号为:%d\n',N_pie); hours_of_N_pie = (N_pie) * (time_interval) - 1 ; %总计经过的小时数 days_of_N_pie = fix(hours_of_N_pie/24); %总天数 real_day_of_N_pie = day + days_of_N_pie; %中间时刻时间不会到四月 real_hour_of_N_pie = hour + rem(hours_of_N_pie , 24); fprintf('中间时刻的世纪时为:'); fprintf('%d年%d月%d日%d时\n',year,month,real_day_of_N_pie,real_hour_of_N_pie); end

        朝朝夕夕,日夜不停。海洋学子在探索海洋的道路上也不会停歇,前有王品先院士,文圣常院士等开创我国海洋学的奠基。未来在维护国家海洋权利,建设海洋强国的道路上,仍然需要大量海洋及其他各个专业学生前赴后继的沉淀。再次希望引起大家的思考和讨论,逐渐完善出可以供大家学习的、标准的算法。

总结

以上是生活随笔为你收集整理的潮汐特征值计算及其预报的全部内容,希望文章能够帮你解决所遇到的问题。

如果觉得生活随笔网站内容还不错,欢迎将生活随笔推荐给好友。