
数值分析
分享数值分析的fortran代码
chder_白南
无话可说
展开
-
matlab生成线性等间隔离散的梯形波
【代码】matlab生成线性等间隔离散的梯形波。原创 2023-12-19 21:59:16 · 967 阅读 · 0 评论 -
fortran实现一维线性插值:与matlab的interp1的线性插值功能一致
program main implicit none integer :: i real*8 :: x(9), y(9), t, z !// x与y的数据来源于matlab的interp1帮助文档, t为被插值节点,见帮助文档的变量xq x = [0.d0,0.785398163397448d0,1.57079632679490d0,2.35619449019235d0,3.14159265358979d0.原创 2021-08-17 23:05:21 · 2373 阅读 · 0 评论 -
matlab寻找一段曲线的极值点
x = [-pi:0.01:pi]';y = sin(10*x);[~,loc1] = findpeaks(y); % 极大值[~,loc2] = findpeaks(max(y)-y); % 极小值原创 2021-04-08 13:22:21 · 9418 阅读 · 0 评论 -
matlab生成任意区间的随机数
使用matlab中的rand函数生成任意区间的随机数公式 r = a + (b-a)*rand(m,n)其中[a,b]是范围,[m,n]是生成的数据形状。比如我想生成[-5,5]之内10个随机数a = -5;b = 5;r = a + (b-a) * rand(10,1); -4.0246 -2.2150 0.4688 4.5751 4.6489 -3.4239 4.7059 4.5717 -0.1462 3.0028.原创 2020-07-26 16:49:20 · 26860 阅读 · 4 评论 -
fortran计算逆矩阵
program main implicit none integer :: n n = 108 call get_invMatrix(n)end program main subroutine get_invMatrix(n) implicit none integer :: n, i, j real*8 :: A(n,n), invA(n,n) ope.原创 2020-07-21 18:42:06 · 7300 阅读 · 1 评论 -
fortran使用MKL函数库求解普通稀疏矩阵与向量的乘积
使用MKL函数库中的mkl_dcsrgemv可以实现普通稀疏矩阵与向量的乘积program main implicit none integer :: n, n1, ja(9), ia(7) real*8 :: a(9), x(6), y(6) character(len=1) :: transa = 'n' !|--------------------------MATRIX--------------.原创 2020-07-18 19:13:13 · 2314 阅读 · 0 评论 -
fortran使用MKL中的GELS求解线性方程组的最小二乘问题
program test_GELS use lapack95 implicit none integer :: i integer, parameter :: m = 4, n = 3 real*8 :: a(4,3), b(4) !|---------------------------------------------------------| .原创 2020-06-04 10:03:07 · 1467 阅读 · 0 评论 -
使用MKL中的有关函数计算共轭梯度法
program ConjugateGradientMethods use blas95 implicit none integer :: i integer, parameter :: n = 5 real(kind=8),allocatable :: A(:,:) real(kind=8),allocatable :: b(:,:) real(kind=8) :: alpha = 0.d0, .原创 2020-05-28 09:23:15 · 697 阅读 · 0 评论 -
自动变步长Simpson积分方法函数
Module mod Implicit none Real(kind=8), parameter :: a = 1.d0, b = 1.02d0 !// 积分区间 Real(kind=8), parameter :: eps = 1.d-15 !// 误差控制 Real(kind=8) :: s !// 最终积分结果 Integer :: n !// 最终分的份数C...原创 2019-12-18 15:24:42 · 1542 阅读 · 0 评论 -
梯形积分与辛普森积分
Program main !// 辛普森积分 Implicit none Real(kind=8), parameter :: a = 0.d0, b = 1.d0 !//积分区间为[a,b], 被积函数为f(x) = x*e^x Integer, parameter :: n = 100 !// 区间等份 Real(kind=8), parameter :: h = (...原创 2019-12-18 15:23:53 · 1274 阅读 · 0 评论 -
龙贝格积分
!// author: luk!// time: 2018/3/20!// purpose: Romberg integral for ln(x) between 1 and 2Module mod Implicit noneContains Subroutine cal_integral( ) Implicit none Integer :: i, j, k ...原创 2019-12-18 15:22:15 · 452 阅读 · 0 评论 -
高斯积分
!// 注: 参考来源http://fcode.cn/algorithm-73-1.html Module Gauss_Legendre !//高斯—勒让德积分高斯点及权重的求解模块 Implicit none Integer, parameter :: n = 11 !// 设置求解高斯点的个数 Integer, parame...原创 2019-12-18 15:21:37 · 758 阅读 · 0 评论 -
复合体型辛普森积分
Program main !// 复合梯形积分和辛普森积分 Implicit none Real(kind=8), parameter :: a = 0.d0, b = 1.d0 !//积分区间为[a,b], 被积函数为f(x) = x*e^x Integer, parameter :: n = 100 !// 区间等份 Real(kind=8), parameter :...原创 2019-12-18 15:21:02 · 417 阅读 · 0 评论 -
自适应辛普森积分
Module Simpson_mod Implicit none Real(kind=8) :: a = 0.d0, b = 1.d0 !// 积分区间 Real(kind=8) :: pi = acos(-1.d0), eps = 1.d-8, s0 Real(kind=8) :: calPi = 0.d0 !// 数值积分解 Contains Real(kind=8...原创 2019-12-18 15:06:59 · 478 阅读 · 0 评论 -
使用Householder方法对矩阵QR分解
Module mod Implicit none Integer, parameter, public :: m = 4, n = 3 !Real(kind=8), public :: A(m,n) = reshape( [1.,2.,2.,-4,-4.,3.,2.,-2,1.,2.,3.,4.],[m,n] ) Real(kind=8), public :: A(m,n) = ...原创 2019-12-18 15:05:44 · 2764 阅读 · 2 评论 -
不完全QR分解
Program QR Implicit none Integer :: i, j Integer, parameter :: m = 4, n = 3 Real(kind=8), parameter :: A(m,n) = reshape( [ 1.,2.,2.,-4.,3.,2.,2.,5.,-2.,6.,-4.,3. ],[m,n] ) Real(kind=8) :: Q...原创 2019-12-18 15:04:46 · 290 阅读 · 0 评论 -
完全QR分解
Program QR_complete Implicit none Integer :: i, j Integer, parameter :: m = 4, n = 2 Real(kind=8) :: A(m,n), Q(m,m), A0(m,n), R(m,n) Real(kind=8) :: t1, t2 a = reshape( [ 1.,2.,2.,-4.,...原创 2019-12-18 15:04:02 · 743 阅读 · 0 评论 -
LM算法求解最小二乘问题
Module VarLM Implicit none !//Nparas为参数个数,Ndata为数据个数,Niters为最大迭代次数 Integer ,Parameter :: Nparas = 2,Ndata = 9,Niters = 50,Fileid = 11 !//临时变量,order=1,继续迭代,否则停止迭代 Integer ...原创 2019-12-18 15:02:41 · 2486 阅读 · 0 评论 -
高斯牛顿法求解最小二乘问题
!// 本代码原理以及测试数据来自wiki: https://en.wikipedia.org/wiki/Gauss%E2%80%93Newton_algorithm#ExampleModule varGS Implicit None !//Nparas为参数个数,ndata为数据个数,Niters为最大迭代次数 Integer, parameter :: Nparas = 2, ...原创 2019-12-16 14:56:48 · 1106 阅读 · 0 评论 -
三次样条插值fortran程序
!// 四类三次样条插值!// 1.自然三次样条插值; 2.钳制三次样条插值; 3.曲率(二阶导数)任意调整的三次样条插值; 4.抛物线端点的三次样条曲线; 5.非纽结三次样条(matlab中的命令为spline)Module CalSpline Implicit none Character(*), parameter :: file1 = 'xy.txt' !// 插值基点文...原创 2019-12-16 14:55:05 · 2568 阅读 · 5 评论 -
牛顿插值fortran程序
!// Newton interpolationModule interpolationContainsSubroutine newtdd( x,y,c ) !// 计算牛顿插值系数 Implicit none Real(kind=8), intent( in ) :: x(:), y(:) Real(kind=8), intent( inout ) :: c(:,:) ...原创 2019-12-16 14:54:28 · 2094 阅读 · 0 评论 -
拉格朗日插值fortran程序
Program LagrangianInterpolation Implicit none Integer :: i, fileid Integer, parameter :: n1 = 51 !// 插值基节点数 Real(kind=8), parameter :: t0 = 0.d0, t1 = 1.d0 !// 插值区间 Real(kind=8) :: x(n1)...原创 2019-12-16 14:53:50 · 3327 阅读 · 1 评论 -
切比雪夫插值fortran程序
!// Chebyshev interpolationModule interpolation ContainsSubroutine newtdd( x,y,c ) !// 计算牛顿插值系数 Implicit none Real(kind=8), intent( in ) :: x(:), y(:) Real(kind=8), intent( inout ) :: c(:...原创 2019-12-16 14:53:07 · 606 阅读 · 0 评论 -
迭代法求解非线性问题
!// 方程组:-u^3 + v = 0!// u^2 + v^2 - 1 = 0!// 即:f1(u,v) = -u^3 + v f2(u,v) = u^2 + v^2 - 1!// 多元牛顿法的求解适用于小型方程组,对于大型不适用,因为需要提前指导导数!// 多变量牛顿法!// x0 = 初始向量!// DF(Xk) * s = -F(Xk)!// Xk+1 ...原创 2019-12-16 14:52:15 · 933 阅读 · 0 评论 -
迭代法求解线性方程组
!// [ 3 -1 0 0 0 0.5 ] [x1] [2.5]!// [ -1 3 -1 0 0.5 0 ] [x2] [1.5]!// [ 0 -1 3 -1 0 0 ] [x3] = [1.0]!// [ 0 0 -1 3 -1 0 ] [x4] [1.0]!// [ 0 0.5 0 -1 3 -1 ] [x5] [1.5]!// [ 0.5 0 0 ...原创 2019-12-16 14:50:47 · 1651 阅读 · 0 评论 -
预条件共轭梯度法求解线性方程组2
!// 在病态矩阵问题上,共轭梯度法比部分主元的高斯消去法还要差!// 这种情况下,可通过预条件得到缓解,主要是将问题转化为良态矩阵系统,再实施CG方法!// 此代码的预条件因子M= ( D + wL ) * Inv(D) * ( D + wU ),为高斯-赛德尔预条件因子!// A = L + D + U,L为A的下三角,U为A的上三角,D为A的对角线元素Program Conjuga...原创 2019-12-16 14:50:01 · 1196 阅读 · 0 评论 -
预条件共轭梯度法求解线性方程组1
!// 在病态矩阵问题上,共轭梯度法比部分主元的高斯消去法还要差!// 这种情况下,可通过预条件得到缓解,主要是将问题转化为良态矩阵系统,再实施CG方法!// 此代码的预条件因子M=D,为雅可比预条件因子,D为A的对角线矩阵Program ConjugateGradientMethods Implicit none Integer :: i, j Integer, parame...原创 2019-12-16 14:49:06 · 2408 阅读 · 3 评论 -
部分主元法求解线性方程组
!//---------------原方程-------------!// x + y + z + w = 10!// 2x + 3y + z + w = 15!// 3x - y + 2z - w = 3!// 4x + y -3z + 2w = 5!//----------------------------------!// 本代码针对高斯消元法进行改...原创 2019-12-16 14:47:50 · 762 阅读 · 0 评论 -
部分主元法求解线性方程组
!//---------------原方程-------------!// x + y + z + w = 10!// 2x + 3y + z + w = 15!// 3x - y + 2z - w = 3!// 4x + y -3z + 2w = 5!//----------------------------------!// 本代码针对高斯消元法进行改...原创 2020-03-20 12:07:43 · 657 阅读 · 0 评论 -
使用置换矩阵求解线性方程组2
!// 说明:至此,第二章有关方程求解的部分就此结束。下节我们正式进入第二章的迭代部分,下期见!!!!// -----------原方程-----------!// 2x + 1y + 5z = 5!// 4x + 4y - 4z = 0!// 1x + 3y + 1z = 6!// ---------------------------!// PA =...原创 2019-11-28 20:18:25 · 602 阅读 · 0 评论 -
使用置换矩阵求解线性方程组1
!// 说明:至此,第二章有关方程求解的部分就此结束。下节我们正式进入第二章的迭代部分,下期见!!!!// -----------原方程-----------!// 2x + 1y + 5z = 5!// 4x + 4y - 4z = 0!// 1x + 3y + 1z = 6!// ---------------------------!// PA =...原创 2019-11-28 20:17:45 · 311 阅读 · 0 评论 -
LU分解求解线性代数方程组
!//---------------原方程-------------!// x + 2y - z = 3!// 2x + y -2z = 3!// -3x + y + z = -6!//----------------------------------!// 本代实现了LU分解,但是并不是所有的矩阵都能进行LU分解!// 后面会介绍PA=LU,可以解决所有问题...原创 2019-11-28 20:16:35 · 1296 阅读 · 0 评论 -
高斯消去法求解线性方程组
!//---------------原方程-------------!// x + y + z + w = 10!// 2x + 3y + z + w = 15!// 3x - y + 2z - w = 3!// 4x + y -3z + 2w = 5!//----------------------------------!// 本代码具有一定的局限性,对...原创 2019-11-28 20:14:33 · 1162 阅读 · 0 评论 -
共轭梯度法求解线性方程组
Program ConjugateGradientMethods Implicit none Integer :: i Integer, parameter :: n = 3 Real(kind=8) :: A(n,n) = [ 4.d0, -2.d0, 2.d0, -2.d0, 2.d0, -4.d0, 2.d0, -4.d0, 11.d0 ] Real(kind=8) :...原创 2019-11-28 20:13:49 · 1941 阅读 · 0 评论 -
楚列斯基分解
!// 首先给出对称正定矩阵的概念!// 如果A'= A,则n*n矩阵A是对称矩阵。如果对于所有向量x/=0,x'Ax > 0,则称矩阵A是正定矩阵!// 对称正定矩阵是非奇异的!// 任意矩阵A的行列式为矩阵对应的特征值的乘积!// 对于对称正定矩阵,可使用楚列斯基分解方法,可分解为: A = R'R,其中R是一个上三角矩阵!// 使用楚列斯基分解,对于对称正定矩阵,和一般的矩...原创 2019-11-28 20:12:23 · 3429 阅读 · 1 评论 -
使用切线法求解方程零点
Program SecantMethod Implicit none Real(kind=8) :: x0 = 0.d0, x1 = 1.d0 Real(kind=8) :: xtmp Real(kind=8), parameter :: eps = 1.d-8 Real(kind=8), external :: func Integer, parameter :: n...原创 2019-11-28 20:12:54 · 958 阅读 · 0 评论 -
使用牛顿迭代法求解方程零点
Program NewtonIterationMethod Implicit none Real(kind=8), external :: func, DerivativeFunc Real(kind=8), parameter :: eps = 1.d-8 Real(kind=8) :: x0, xtmp Integer, parameter :: nloop = 100...原创 2019-11-28 20:09:48 · 2174 阅读 · 0 评论 -
使用不动点迭代方法求解方程的零点
Program FixedPointIteration !// 使用不动点迭代方法求解方程 x^3 + x - 1 = 0 的零点 Implicit none Real(kind=8), parameter :: eps = 1.d-8 Real(kind=8), external :: func Real(kind=8) :: xtemp, x0 = 0.5d0 !...原创 2019-11-28 20:09:00 · 2020 阅读 · 0 评论 -
二分法求解方程的零点
Program BisectionMethod !// 二分法求解方程f(t) = cos(t) - t在区间[0,1]上的根 Implicit none Real(kind=8), parameter :: eps = 1d-12 !// 精度控制 Real(kind=8) :: a = 0.d0, b = 1.d0, c, t Real(kind=8) :: fa, fb,...原创 2019-11-28 20:08:08 · 523 阅读 · 0 评论