Ruby

数值计算编程作业(Ruby版)

Posted on

内容包括

  1. 二分法求解二元一次方程的根
  2. 牛顿法求解一元多次方程的根
  3. 顺序消去法解线性方程组
  4. 列选主元消去法解线性方程组
  5. 全选主元消去法解线性方程组
  6. Doolittle分解线性方程组
  7. Crout分解线性方程组
  8. 平方根法求线性方程组
  9. 拉格朗日插值法
  10. 牛顿插值法
  11. 最小二乘法——线性拟合
  12. 变步长梯形求积算法计算积分
  13. 龙贝格算法计算积分
  14. 改进欧拉算法求常微分方程的数值解
  15. 四阶龙格-库塔法求常微分方程的数值解

累!!!!!!!不详细介绍了

#!/usr/bin/env ruby
# encoding: utf-8

$eps=1e-8
$inf=0x3f3f3f3f

class TwoSplit
  def initialize
    @a,@b,@c,@d=2,4,2,18
  end

  def getVal(x)
    @a*x*x+@b*x+@c-@d
  end
  private :getVal

  def three_split
    le,ri,mid=-$inf,$inf,1
    while le+$eps<ri
      mid=(le+ri)/2
      if(getVal(mid)+$eps<0)
        le=mid+$eps
      else
        ri=mid-$eps
      end
    end
    mid
  end
  private :three_split

  def two_split(le,ri)
    mid=1
    while le+$eps<ri
      mid=(le+ri)/2
      if(getVal(mid)<0)
        le=mid+$eps
      else
        ri=mid-$eps
      end
    end
    mid
  end
  private :two_split

  def calculate_ans
    @a,@b,@c,@d=-@a,-@b,-@c,-@d if (@a<0)
    k=three_split
    if(getVal(k)>0)
      puts 'No Answer'
    else
      rans=two_split(k,1.0*$inf)
      lans=2*k-rans
      puts '测试方程为2x^2+4x+2=18'
      printf("方程解为:%.2f\n",lans)
      printf(" %.2f",rans) if((lans-rans).abs>$eps*10)
      puts ''
    end
  end
end

class NewtonIterate
  def initialize
    @a,@b,@c,@d=-2,-5,1,3
  end

  def get_val(x)
    @a*x*x+@b*x+@c-@d
  end

  def get_der(x)
    2*@a*x+@b
  end

  def newton(val)
    pre=val
    val=pre-get_val(pre)/get_der(pre)
    while((pre-val).abs>$eps)
      pre=val
      val=pre-get_val(pre)/get_der(pre)
    end
    val
  end

  def calculate_ans
    if(@b*@b-4*@a*@c<0)
      puts 'No Answer'
    else
      lans=newton(-1.0*$inf)
      rans=newton(1.0*$inf)
      puts '求解方程为-2x^2-5x+1=3'
      printf("方程解为:%.2f",lans);
      printf(" %.2f",rans) if((lans-rans).abs>$eps)
      puts ''
    end
  end
end

class ResolveMat
  def initialize
    @mat=6.times.map{ [0] * 6 }
    @l=6.times.map{ [0] * 6 }
    @u=6.times.map{ [0] * 6 }
    @x,@y=[0]*6,[0]*6
    @b=[1,2,3,4]
    @mat[0]=[2,10,4,-3]
    @mat[1]=[-3,-4,-12,13]
    @mat[2]=[1,2,6,-4]
    @mat[3]=[4,14,9,-13]
    @n=4
  end

  def print_mat(mat)
    yield
    @n.times do |i|
      print '['
      @n.times do |j|
        if((mat[i][j]-0).abs > $eps)
          printf("%8.3f ", mat[i][j]);
        else
          printf("         ");
        end
      end
      puts ']'
    end
    puts ''
  end

  def doolittle
    @n.times { |i| @l[i][i]=1 }

    @n.times do |k|
      k.upto(@n-1) do |j|
        @u[k][j]=@mat[k][j]
        k.times { |i| @u[k][j]-=(@l[k][i]*@u[i][j]) }
      end

      (k+1).upto(@n-1) do |i|
        @l[i][k]=@mat[i][k]
        k.times { |j| @l[i][k]-=(@l[i][j]*@u[j][k]) }
        @l[i][k]/=1.0*@u[k][k]
      end
    end

    @n.times do |i|
      @y[i]=@b[i]
      i.times { |j| @y[i]-=(@l[i][j]*@y[j]) }
    end

    (@n-1).downto(0) do |i|
      @x[i]=@y[i]
      (i+1).upto(@n-1) { |j| @x[i]-=(@u[i][j]*@x[j]) }
      @x[i]/=1.0*@u[i][i]
    end
    puts 'Doolittle分解矩阵'
    print_mat(@mat) {puts '原矩阵为:'}
    print_mat(@l) {puts '矩阵L为:'}
    print_mat(@u) {puts '矩阵U为:'}
  end

  def crout
    @n.times { |i| @u[i][i]=1 }

    @n.times do |k|
      k.upto(@n-1) do |i|
        @l[i][k]=@mat[i][k]
        k.times { |j| @l[i][k]-=(@l[i][j]*@u[j][k]) }
      end
      (k+1).upto(@n-1) do |i|
        @u[k][i]=@mat[k][i]
        k.times { |j| @u[k][i]-=(@l[k][j]*@u[j][i]) }
        @u[k][i]/=1.0*@l[k][k]
      end
    end

    @n.times do |i|
      @y[i]=@b[i]
      i.times { |j| @y[i]-=(@l[i][j]*@y[j]) }
      @y[i]/=1.0*@l[i][i]
    end
    (@n-1).downto(0) do |i|
      @x[i]=@y[i]
      (i+1).upto(@n-1) { |j| @x[i]-=(@u[i][j]*@x[j]) }
      # @x[i]/=@u[i][i]
    end
    puts 'Crout分解矩阵'
    print_mat(@mat) {puts '原矩阵为:'}
    print_mat(@l) {puts '矩阵L为:'}
    print_mat(@u) {puts '矩阵U为:'}
  end

  def square #待验证
    t=6.times.map { [0]*6 }
    lt=6.times.map { [0]*6 }
    d=6.times.map { [0]*6 }

    @n.times { |i| lt[i][i]=1 }
    t[0][0]=@mat[0][0]

    1.upto(@n-1) do |i|
      t[i][0] = @mat[i][0]
      lt[0][i] = @mat[0][i] / t[0][0]
    end

    1.upto(@n-1) do |i|
      1.upto(i) do |j|
        t[i][0]=@mat[i][0]
        j.times { |k| t[i][j]-=t[i][k]*lt[k][j] }
      end

      (i+1).upto(@n-1) do |j|
        lt[i][j] = @mat[i][j]
        i.times {|k| lt[i][j] = lt[i][j] - t[i][k] * lt[k][j] }
        lt[i][j] = lt[i][j] / t[i][i]
      end
    end

    @n.times do |i|
       @n.times do |j|
         @l[i][j]=lt[i][j]
         d[i][j]=0 if i!=j
         d[i][j]=t[i][j] if i==j
       end
    end

    @y[0]=@b[0]
    1.upto(@n-1) do |i|
      i.times.reduce(0.0) { |sum,k| @y[i]= @b[i]-(sum+=lt[k][i]*@y[k]) }
      # m=0.0
    end

    @x[@n-1]=@y[@n-1]/t[@n-1][@n-1]
    (@n-1).downto(0) do |i|
      (i+1).upto(@n-1).reduce(0.0)  { |sum,k|  @x[i]=@y[i]/d[i][i] - (sum+=lt[i][k]*@x[k]) }
    end

    print_mat(@mat) { puts '所要求解的线性方程为:'}
    @n.times { |i| printf("x[%d]=%.3f \n", i+1, @x[i]); }
  end
end

class Order

  def initialize
    @mat = Array.new(4,[])
    @mat[1]=[0,1,2,3,1]
    @mat[2]=[0,5,4,10,0]
    @mat[3]=[0,3,-0.1,1,2]
    @n=3
  end

  def print_mat
    yield
    @n.times do | i |
      (@n+1).times do | j |
        printf("%.2f ",@mat[i+1][j+1])
      end
      puts ''
    end
  end
  private :print_mat

  def print_ans
    print "解为:"
    @n.times do | k |
      printf("%.2f ",@mat[k+1][@n+1])
    end
    printf "\n\n"
  end
  private :print_ans

  def calculate_answer
    @mat[@n][@n+1] /= @mat[@n][@n]
    (@n-1).downto(1) do |i|
      @mat[i][@n+1]-=(i+1..@n).reduce(0) {|sum,j| sum+=@mat[i][j]*@mat[j][@n+1]}
    end
  end
  private :calculate_answer

  def calculate_mat(i)
    (i+1).upto(@n+1) do | j |
      @mat[i][j]/=(1.0*@mat[i][i])
    end
    @mat[i][i]=1
    (i+1).upto(@n) do | k |
      (i+1).upto(@n+1) do | j |
        @mat[k][j]-=@mat[k][i]*@mat[i][j]
      end
    end
  end
  private :calculate_mat

  def original_order
    print_mat {puts "顺序消去法:\n原矩阵为:"}
    1.upto(@n-1) do | i |
      calculate_mat(i)
    end
    print_mat {puts '所得上三角矩阵为:'}
    calculate_answer
    print_ans
  end

  def col_order
    print_mat {puts "列主元素消去法:\n原矩阵为:"}
    1.upto(@n-1) do | i |
      now,maxx=i,@mat[i][i]
      (i+1).upto(@n) do | j |
        if(@mat[j][i]>maxx)
          maxx,now=@mat[j][i],j
        end
      end

      if i!=now
        (i..@n+1).each do | k |
         @mat[i][k],@mat[now][k] = @mat[now][k],@mat[i][k]
        end
      end
    calculate_mat(i)
    end
    print_mat {puts '所得上三角矩阵为:'}
    calculate_answer
    print_ans
  end

  def all_order
    print_mat {puts "全主元素消去法:\n原矩阵为:"}
    1.upto(@n-1) do | i |
      nx,ny,maxx=i,i,@mat[i][i]

      i.upto(@n) do | j |
        i.upto(@n) do | t |
          if(@mat[j][t]>maxx)
            maxx,nx,ny=@mat[j][t],j,t
          end
        end
      end

      if i!=nx #row exchange
        i.upto(@n+1) do | k |
          @mat[i][k],@mat[nx][k] = @mat[nx][k],@mat[i][k]
        end
      end

      if i!=ny # col exchange
        i.upto(@n) do | k |
          @mat[k][i],@mat[k][ny] = @mat[k][ny],@mat[k][i]
        end
      end

      calculate_mat(i)
    end
    print_mat {puts '所得上三角矩阵为:'}
    calculate_answer
    print_ans
  end

end

class VarStep
  def initialize
    @val = (1 + getVal(1))/2
    @pre,@n = 0x3f3f3f3f,1
  end

  def getVal(x)
    Math.sin(x) / x
  end

  def func(n,pre)
    sum,x = 0.0,0
    n.times do
      x = x + 1.0/(2*n)
      sum = sum + getVal(x)
      x = x + 1.0/(2*n)
    end
    @pre/2.0 + sum/(2*n)
  end

  def calculate_ans
    while (@val-@pre).abs > $eps
      @pre = @val
      @val = func(@n,@val)
      @n*=2
    end
    puts '积分函数为sin(x)/x,积分区间为[0,1]'
    printf("积分结果为%.8f\n",@val)
  end

end

class FourthRomberg
  def initialize
    @a,@b,@n=1,5,5
    @x,@y,@k=[],[],[]
    @h=1.0*(@b-@a)/@n
    @x[0],@y[0]=@a,17
  end

  def der_func(x,y)
    x-y+1
  end
  private :der_func

  def calculate_ans
    1.upto(@n) do |i|
      @x[i-1]=@x[0]+(i-1)*@h
      h1=@h/2.0
      @k[1]=der_func(@x[i-1],@y[i-1])
      @k[2]=der_func(@x[i-1]+h1,@y[i-1]+h1*@k[1])
      @k[3]=der_func(@x[i-1]+h1,@y[i-1]+h1*@k[2])
      @k[4]=der_func(@x[i-1]+@h,@y[i-1]+@h*@k[3])
      @y[i]=@y[i-1]+@h*(@k[1]+2*@k[2]+2*@k[3]+@k[4])/6.0
    end

    puts "a=1,b=5\nn=5\ny0=17\nf(x)=x-y+1\n"
    @n.times do |i|
      printf("y[%d]: %.2f\n",i+1,@y[i+1])
    end
  end
end

class ImproveEuler
  def initialize
    @up,@down,@n = 0.5,0,5
    @h = (@up-@down) / @n
    @xn = @down
    @yn = 1
  end

  def der_func(x,y)
    x - y + 1
  end

  def calculate_ans
    print "f(x,y)=x-y+1\n积分区间为[0,0.5]\ny0=1\n区间5等分\n"
    @n.times do | i |
      xl = @xn + @h
      ylb = @yn + @h * der_func(@xn,@yn)
      yl = @yn + @h / 2 * (der_func(@xn,@yn)+der_func(xl,ylb))
      printf("x%d=%f,y%d=%f\n", i, xl, i, yl)
      @xn,@yn = xl,yl
    end
  end
end

class Lagrance
  def initialize
    @n,@m,@ans=2,195,0
    @x=[121,144,100]
    @y=[11,12,10]
  end

  def calculate_ans
    0.upto(@n) do |i|
      key=1.0*@y[i]
      0.upto(@n) do |j|
        next if(i==j)
        key*=(1.0*(@m-@x[j])/(@x[i]-@x[j]))
      end
      @ans+=key
    end
    puts '给定3个坐标:'
    puts 'x[0]:121,y[0]:11'
    puts 'x[1]:144,y[1]:12'
    puts 'x[2]:100,y[2]:10'
    printf("询问195的y值为%.2f\n",@ans)
  end
end

class NewtonInterpolate
  def initialize
    @n=4
    @form=[0.0]*6
    @x=[10,9,3,11,2]
    @y=[100,80,9,123,4]
  end

  def calculate_ans
    @n.times { |i| @form[i]=@y[i] }
    (@n+1).times do |i|
      @n.downto(i+1) { |j| @form[j] = (@form[j]-@form[j - 1]) / (@x[j] - @x[j - 1 - i]) }
    end
    tmp,hx=1,6.42
    newton=@form[0]
    @n.times do |i|
      tmp*=(hx-@x[i])
      newton+=tmp*@form[i+1]
    end
    @n.times { |i| printf("第%d个点%.2f %.2f\n",i+1,@x[i],@y[i]) }
    printf("F(%.2f) = %.2f\n", hx, newton)
  end
end

class LeastSquare
  def initialize
    @x=[10,9,5,6,2]
    @y=[99,81,24,37,4]
    @a,@b=0.0,0.0
    @n=5
  end

  def get_val(x)
    @a*x+@b
  end

  def calculate_ans
    t=[0.0]*4
    @n.times do |i|
      t[0]+=@x[i]*@x[i]
      t[1]+=@x[i]
      t[2]+=@x[i]*@y[i]
      t[3]+=@y[i]
    end
    @a=(t[2]*@n-t[1]*t[3])/(t[0]*@n-t[1]*t[1])
    @b=(t[0]*t[3]-t[1]*t[2])/(t[0]*@n-t[1]*t[1])
    printf("所得函数为%dx%d\n",@a,@b)
    puts '测试数:12'
    puts get_val(12)
  end
end

class Romberg
  @@pi=4.0*Math.atan(1.0)
  def initialize
    @a,@b=0.000001,1.0
    @n,@h=20,Float(@b-@a)
    @t=20.times.map{ [0] * 2 }
    @t[0][1]=trapezium
    @n*=2
  end

  def get_val(x)
    Math.sin(x)/x
  end
  private :get_val

  def trapezium
    @h=(@b-@a)/@n
    sum = 1.upto(@n-1).reduce(0.0) {|tmp,i| tmp+=get_val(@a+i*@h) }
    sum+=(get_val(@a)+get_val(@b))/2.0
    @h*sum
  end

  def calculate_ans

    puts 'sin(x)/x在[0,1]区间内的积分测试:'
    (1..9).each do |m|
      m.times { |i| @t[i][0]=@t[i][1] }
      @t[0][1]=trapezium
      @n*=2
      m.times do |i|
        @t[i+1][1]=@t[i][1]+(@t[i][1]-@t[i][0]) / ((4**m)-1)
      end
      if (@t[m-1][1]-@t[m][1]).abs < $eps
        printf("计算的数为:%.6f\n", @t[m][1])
        return @t[m][1]
      end
    end
    puts 'No Answer'
  end
end

def nodekey
  while true
    puts '<===================================>'
    puts '1.二分法求解二元一次方程的根'
    puts '2.牛顿法求解一元多次方程的根'
    puts '3.顺序去法解线性方程组'
    puts '4.列选主元消去法解线性方程组'
    puts '5.全选主元消去法解线性方程组'
    puts '6.Doolittle分解线性方程组'
    puts '7.Crout分解线性方程组'
    puts '8.平方根法求线性方程组'
    puts '9.拉格朗日插值法'
    puts '10.牛顿插值法'
    puts '11.最小二乘法——线性拟合'
    puts '12.变步长梯形求积算法计算积分'
    puts '13.龙贝格算法计算积分'
    puts '14.改进欧拉算法求常微分方程的数值解'
    puts '15.四阶龙格-库塔法求常微分方程的数值解'
    puts ' 0.退出'
    puts '<====================================>'
    print "请选择你要选择的操作:"
    op=gets.to_i
    case op
      when 0
        break
      when 1
        TwoSplit.new.calculate_ans
      when 2
        NewtonIterate.new.calculate_ans
      when 3
        Order.new.original_order
      when 4
        Order.new.col_order
      when 5
        Order.new.all_order
      when 6
        ResolveMat.new.doolittle
      when 7
        ResolveMat.new.crout
      when 8
        ResolveMat.new.square
      when 9
        Lagrance.new.calculate_ans
      when 10
        NewtonInterpolate.new.calculate_ans
      when 11
        LeastSquare.new.calculate_ans
      when 12
        VarStep.new.calculate_ans
      when 13
        Romberg.new.calculate_ans
      when 14
        ImproveEuler.new.calculate_ans
      when 15
        FourthRomberg.new.calculate_ans
      else
        puts '请输入0-15的选项,选择0为退出'
    end
  end
end

nodekey
哈夫曼树

哈夫曼编码

Posted on

传送门

OJ小作业 很久以前第一次听说哈夫曼树问了一下煞笔室友 听他很不屑地说 只是一个贪心而已 我也就没怎么在意

最近李总上课点名不得不去上一下 可惜也只是到那里发了一下呆 点了一下到就过去了 他在讲哈夫曼编码的时候我也没去听

搞得我刚开始看到这题的时候有点小蒙蔽 还好数据比较水 不然估计要卡很长时间 (因为完全是我根据以前的知识xjb搞的 还好还好是一A

下面是源码+注解 (因为完全是自己看了远离后搞得 可能有点烦........

#include <algorithm>
#include <cmath>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <queue>
#include <stack>

#define root tree[id]
#define leson tree[ldir]
#define rison tree[rdir]

using namespace std;
const int maxn = 104;
int len, all_time;
int tree_root;//根节点
char to_c[maxn];//由节点的值转化为字符
char tar[maxn * 100];
char tar_ch;//当前目标字符
char output[maxn];//输出数组
bool flag;//标记  如果已经输出就直接返回

struct node {
    int lson, rson;
    int tim;//结点生成时间  所有原结点初始化为0
    int order;//节点编号
    int val;//结点的值
    bool operator<(const node& a) const
    {
        if (val != a.val)
            return val > a.val;
        if (tim != a.tim)
            return tim > a.tim;
        return order > a.order;
    }//按照题意出队
};

priority_queue<node> que;
node tree[maxn << 2];

int CreatTree()
{
    node ta, tb;
    int id;
    while (que.size() > 1) {
        ta = que.top();
        que.pop();
        tb = que.top();
        que.pop();
        id = len + all_time;
        root.lson = ta.order;
        root.rson = tb.order;
        root.val = ta.val + tb.val;
        root.order = id;
        root.tim = all_time++;
        que.push(root);
        // printf("%d %d\n",ta.order,tb.order);
    }//取出两个节点并合并再入队   最后剩余的便是根节点
    ta = que.top();
    que.pop();
    return ta.order;
}

void lock(int id, int dep)
{
    if (flag)//若已经输出过了   就直接返回
        return;
    int ldir = root.lson, rdir = root.rson;
    if (ldir) {
        output[dep] = '0';
        lock(ldir, dep + 1);
    }
    if (rdir) {
        output[dep] = '1';
        lock(rdir, dep + 1);
    }
    if (!ldir && !rdir && to_c[root.order] == tar_ch) {  //确保该结点为叶子结点也可以省点时间
        for (int i = 0; i < dep; i++)
            putchar(output[i]);
        flag = true;//标记为true
    }
}

int unlock(int index)
{
    int id = tree_root;
    while (root.lson || root.rson) {
        if (tar[index] == '0')  // 0 向左 1向右
            id = root.lson;
        else
            id = root.rson;
        index++;
    }
    putchar(to_c[root.order]);
    return index;//返回转化后的index
}

int main()
{
    while (scanf("%s", to_c + 1) != EOF) {
        len = strlen(to_c + 1);
        for (int id = 1; id <= len; id++) {
            scanf("%d", &root.val);
            root.order = id;
            root.tim = 0;
            root.lson = root.rson = 0;
            que.push(root);
        }
        all_time = 1;
        tree_root = CreatTree();
        // printf("root: %d\n",tree_root);
        // for(int id=1;id<=tree_root;id++)
        //    printf("order: %d  val: %d  lson: %d  rson:%d\n",root.order,root.val,root.lson,root.rson); //可取消注释查看建树是否正确
        scanf("%s", tar);
        int tar_len = strlen(tar);
        for (int i = 0; i < tar_len; i++) {
            tar_ch = tar[i];
            flag = false;
            lock(tree_root, 0);
        }
        puts("");
        scanf("%s", tar);
        tar_len = strlen(tar);
        int i = 0;
        while (i < tar_len)
            i = unlock(i);
        puts("");
    }
    return 0;
}