Matlab一度被认为是最专业的数值计算工具之一,相信许多同学都或多或少用过这个工具。相比而言,Python作为一种胶水式的语言,其设计之初就不是为科学计算服务的。之前也看到许多人在吐槽说用Python去复现一些计算过程时经常失败,因此(包括本人)也怀疑过是Python本身数值精度不够导致的。那么Python的精度究竟如何,本文就来一探究竟。为了方便,我们就用线性方程组的求解来对比这一事实。

1、实验设计
  • 基本思路:

本文考虑的线性方程组:
Ax=bAx=b Ax=b

#mermaid-svg-KlQCVM1aHWdwKg9a .label {font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family);fill: #333;color: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .label text {fill: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .node rect, #mermaid-svg-KlQCVM1aHWdwKg9a .node circle, #mermaid-svg-KlQCVM1aHWdwKg9a .node ellipse, #mermaid-svg-KlQCVM1aHWdwKg9a .node polygon, #mermaid-svg-KlQCVM1aHWdwKg9a .node path {fill: #ECECFF;stroke: #9370DB;stroke-width: 1px; }#mermaid-svg-KlQCVM1aHWdwKg9a .node .label {text-align: center;fill: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .node.clickable {cursor: pointer; }#mermaid-svg-KlQCVM1aHWdwKg9a .arrowheadPath {fill: #333333; }#mermaid-svg-KlQCVM1aHWdwKg9a .edgePath .path {stroke: #333333;stroke-width: 1.5px; }#mermaid-svg-KlQCVM1aHWdwKg9a .flowchart-link {stroke: #333333;fill: none; }#mermaid-svg-KlQCVM1aHWdwKg9a .edgeLabel {background-color: #e8e8e8;text-align: center; } #mermaid-svg-KlQCVM1aHWdwKg9a .edgeLabel rect {opacity: 0.9; } #mermaid-svg-KlQCVM1aHWdwKg9a .edgeLabel span {color: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .cluster rect {fill: #ffffde;stroke: #aaaa33;stroke-width: 1px; }#mermaid-svg-KlQCVM1aHWdwKg9a .cluster text {fill: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a div.mermaidTooltip {position: absolute;text-align: center;max-width: 200px;padding: 2px;font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family);font-size: 12px;background: #ffffde;border: 1px solid #aaaa33;border-radius: 2px;pointer-events: none;z-index: 100; }#mermaid-svg-KlQCVM1aHWdwKg9a .actor {stroke: #CCCCFF;fill: #ECECFF; }#mermaid-svg-KlQCVM1aHWdwKg9a text.actor > tspan {fill: black;stroke: none; }#mermaid-svg-KlQCVM1aHWdwKg9a .actor-line {stroke: grey; }#mermaid-svg-KlQCVM1aHWdwKg9a .messageLine0 {stroke-width: 1.5;stroke-dasharray: none;stroke: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .messageLine1 {stroke-width: 1.5;stroke-dasharray: 2, 2;stroke: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a #arrowhead path {fill: #333;stroke: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .sequenceNumber {fill: white; }#mermaid-svg-KlQCVM1aHWdwKg9a #sequencenumber {fill: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a #crosshead path {fill: #333;stroke: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .messageText {fill: #333;stroke: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .labelBox {stroke: #CCCCFF;fill: #ECECFF; }#mermaid-svg-KlQCVM1aHWdwKg9a .labelText,#mermaid-svg-KlQCVM1aHWdwKg9a .labelText > tspan {fill: black;stroke: none; }#mermaid-svg-KlQCVM1aHWdwKg9a .loopText,#mermaid-svg-KlQCVM1aHWdwKg9a .loopText > tspan {fill: black;stroke: none; }#mermaid-svg-KlQCVM1aHWdwKg9a .loopLine {stroke-width: 2px;stroke-dasharray: 2, 2;stroke: #CCCCFF;fill: #CCCCFF; }#mermaid-svg-KlQCVM1aHWdwKg9a .note {stroke: #aaaa33;fill: #fff5ad; }#mermaid-svg-KlQCVM1aHWdwKg9a .noteText,#mermaid-svg-KlQCVM1aHWdwKg9a .noteText > tspan {fill: black;stroke: none; }#mermaid-svg-KlQCVM1aHWdwKg9a .activation0 {fill: #f4f4f4;stroke: #666; }#mermaid-svg-KlQCVM1aHWdwKg9a .activation1 {fill: #f4f4f4;stroke: #666; }#mermaid-svg-KlQCVM1aHWdwKg9a .activation2 {fill: #f4f4f4;stroke: #666; }#mermaid-svg-KlQCVM1aHWdwKg9a .mermaid-main-font {font-family: "trebuchet ms", verdana, arial;font-family: var(--mermaid-font-family); }#mermaid-svg-KlQCVM1aHWdwKg9a .section {stroke: none;opacity: 0.2; }#mermaid-svg-KlQCVM1aHWdwKg9a .section0 {fill: rgba(102, 102, 255, 0.49); }#mermaid-svg-KlQCVM1aHWdwKg9a .section2 {fill: #fff400; }#mermaid-svg-KlQCVM1aHWdwKg9a .section1, #mermaid-svg-KlQCVM1aHWdwKg9a .section3 {fill: white;opacity: 0.2; }#mermaid-svg-KlQCVM1aHWdwKg9a .sectionTitle0 {fill: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .sectionTitle1 {fill: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .sectionTitle2 {fill: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .sectionTitle3 {fill: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .sectionTitle {text-anchor: start;font-size: 11px;text-height: 14px;font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family); }#mermaid-svg-KlQCVM1aHWdwKg9a .grid .tick {stroke: lightgrey;opacity: 0.8;shape-rendering: crispEdges; } #mermaid-svg-KlQCVM1aHWdwKg9a .grid .tick text {font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family); }#mermaid-svg-KlQCVM1aHWdwKg9a .grid path {stroke-width: 0; }#mermaid-svg-KlQCVM1aHWdwKg9a .today {fill: none;stroke: red;stroke-width: 2px; }#mermaid-svg-KlQCVM1aHWdwKg9a .task {stroke-width: 2; }#mermaid-svg-KlQCVM1aHWdwKg9a .taskText {text-anchor: middle;font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family); }#mermaid-svg-KlQCVM1aHWdwKg9a .taskText:not([font-size]) {font-size: 11px; }#mermaid-svg-KlQCVM1aHWdwKg9a .taskTextOutsideRight {fill: black;text-anchor: start;font-size: 11px;font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family); }#mermaid-svg-KlQCVM1aHWdwKg9a .taskTextOutsideLeft {fill: black;text-anchor: end;font-size: 11px; }#mermaid-svg-KlQCVM1aHWdwKg9a .task.clickable {cursor: pointer; }#mermaid-svg-KlQCVM1aHWdwKg9a .taskText.clickable {cursor: pointer;fill: #003163 !important;font-weight: bold; }#mermaid-svg-KlQCVM1aHWdwKg9a .taskTextOutsideLeft.clickable {cursor: pointer;fill: #003163 !important;font-weight: bold; }#mermaid-svg-KlQCVM1aHWdwKg9a .taskTextOutsideRight.clickable {cursor: pointer;fill: #003163 !important;font-weight: bold; }#mermaid-svg-KlQCVM1aHWdwKg9a .taskText0, #mermaid-svg-KlQCVM1aHWdwKg9a .taskText1, #mermaid-svg-KlQCVM1aHWdwKg9a .taskText2, #mermaid-svg-KlQCVM1aHWdwKg9a .taskText3 {fill: white; }#mermaid-svg-KlQCVM1aHWdwKg9a .task0, #mermaid-svg-KlQCVM1aHWdwKg9a .task1, #mermaid-svg-KlQCVM1aHWdwKg9a .task2, #mermaid-svg-KlQCVM1aHWdwKg9a .task3 {fill: #8a90dd;stroke: #534fbc; }#mermaid-svg-KlQCVM1aHWdwKg9a .taskTextOutside0, #mermaid-svg-KlQCVM1aHWdwKg9a .taskTextOutside2 {fill: black; }#mermaid-svg-KlQCVM1aHWdwKg9a .taskTextOutside1, #mermaid-svg-KlQCVM1aHWdwKg9a .taskTextOutside3 {fill: black; }#mermaid-svg-KlQCVM1aHWdwKg9a .active0, #mermaid-svg-KlQCVM1aHWdwKg9a .active1, #mermaid-svg-KlQCVM1aHWdwKg9a .active2, #mermaid-svg-KlQCVM1aHWdwKg9a .active3 {fill: #bfc7ff;stroke: #534fbc; }#mermaid-svg-KlQCVM1aHWdwKg9a .activeText0, #mermaid-svg-KlQCVM1aHWdwKg9a .activeText1, #mermaid-svg-KlQCVM1aHWdwKg9a .activeText2, #mermaid-svg-KlQCVM1aHWdwKg9a .activeText3 {fill: black !important; }#mermaid-svg-KlQCVM1aHWdwKg9a .done0, #mermaid-svg-KlQCVM1aHWdwKg9a .done1, #mermaid-svg-KlQCVM1aHWdwKg9a .done2, #mermaid-svg-KlQCVM1aHWdwKg9a .done3 {stroke: grey;fill: lightgrey;stroke-width: 2; }#mermaid-svg-KlQCVM1aHWdwKg9a .doneText0, #mermaid-svg-KlQCVM1aHWdwKg9a .doneText1, #mermaid-svg-KlQCVM1aHWdwKg9a .doneText2, #mermaid-svg-KlQCVM1aHWdwKg9a .doneText3 {fill: black !important; }#mermaid-svg-KlQCVM1aHWdwKg9a .crit0, #mermaid-svg-KlQCVM1aHWdwKg9a .crit1, #mermaid-svg-KlQCVM1aHWdwKg9a .crit2, #mermaid-svg-KlQCVM1aHWdwKg9a .crit3 {stroke: #ff8888;fill: red;stroke-width: 2; }#mermaid-svg-KlQCVM1aHWdwKg9a .activeCrit0, #mermaid-svg-KlQCVM1aHWdwKg9a .activeCrit1, #mermaid-svg-KlQCVM1aHWdwKg9a .activeCrit2, #mermaid-svg-KlQCVM1aHWdwKg9a .activeCrit3 {stroke: #ff8888;fill: #bfc7ff;stroke-width: 2; }#mermaid-svg-KlQCVM1aHWdwKg9a .doneCrit0, #mermaid-svg-KlQCVM1aHWdwKg9a .doneCrit1, #mermaid-svg-KlQCVM1aHWdwKg9a .doneCrit2, #mermaid-svg-KlQCVM1aHWdwKg9a .doneCrit3 {stroke: #ff8888;fill: lightgrey;stroke-width: 2;cursor: pointer;shape-rendering: crispEdges; }#mermaid-svg-KlQCVM1aHWdwKg9a .milestone {transform: rotate(45deg) scale(0.8, 0.8); }#mermaid-svg-KlQCVM1aHWdwKg9a .milestoneText {font-style: italic; }#mermaid-svg-KlQCVM1aHWdwKg9a .doneCritText0, #mermaid-svg-KlQCVM1aHWdwKg9a .doneCritText1, #mermaid-svg-KlQCVM1aHWdwKg9a .doneCritText2, #mermaid-svg-KlQCVM1aHWdwKg9a .doneCritText3 {fill: black !important; }#mermaid-svg-KlQCVM1aHWdwKg9a .activeCritText0, #mermaid-svg-KlQCVM1aHWdwKg9a .activeCritText1, #mermaid-svg-KlQCVM1aHWdwKg9a .activeCritText2, #mermaid-svg-KlQCVM1aHWdwKg9a .activeCritText3 {fill: black !important; }#mermaid-svg-KlQCVM1aHWdwKg9a .titleText {text-anchor: middle;font-size: 18px;fill: black;font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family); }#mermaid-svg-KlQCVM1aHWdwKg9a g.classGroup text {fill: #9370DB;stroke: none;font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family);font-size: 10px; } #mermaid-svg-KlQCVM1aHWdwKg9a g.classGroup text .title {font-weight: bolder; }#mermaid-svg-KlQCVM1aHWdwKg9a g.clickable {cursor: pointer; }#mermaid-svg-KlQCVM1aHWdwKg9a g.classGroup rect {fill: #ECECFF;stroke: #9370DB; }#mermaid-svg-KlQCVM1aHWdwKg9a g.classGroup line {stroke: #9370DB;stroke-width: 1; }#mermaid-svg-KlQCVM1aHWdwKg9a .classLabel .box {stroke: none;stroke-width: 0;fill: #ECECFF;opacity: 0.5; }#mermaid-svg-KlQCVM1aHWdwKg9a .classLabel .label {fill: #9370DB;font-size: 10px; }#mermaid-svg-KlQCVM1aHWdwKg9a .relation {stroke: #9370DB;stroke-width: 1;fill: none; }#mermaid-svg-KlQCVM1aHWdwKg9a .dashed-line {stroke-dasharray: 3; }#mermaid-svg-KlQCVM1aHWdwKg9a #compositionStart {fill: #9370DB;stroke: #9370DB;stroke-width: 1; }#mermaid-svg-KlQCVM1aHWdwKg9a #compositionEnd {fill: #9370DB;stroke: #9370DB;stroke-width: 1; }#mermaid-svg-KlQCVM1aHWdwKg9a #aggregationStart {fill: #ECECFF;stroke: #9370DB;stroke-width: 1; }#mermaid-svg-KlQCVM1aHWdwKg9a #aggregationEnd {fill: #ECECFF;stroke: #9370DB;stroke-width: 1; }#mermaid-svg-KlQCVM1aHWdwKg9a #dependencyStart {fill: #9370DB;stroke: #9370DB;stroke-width: 1; }#mermaid-svg-KlQCVM1aHWdwKg9a #dependencyEnd {fill: #9370DB;stroke: #9370DB;stroke-width: 1; }#mermaid-svg-KlQCVM1aHWdwKg9a #extensionStart {fill: #9370DB;stroke: #9370DB;stroke-width: 1; }#mermaid-svg-KlQCVM1aHWdwKg9a #extensionEnd {fill: #9370DB;stroke: #9370DB;stroke-width: 1; }#mermaid-svg-KlQCVM1aHWdwKg9a .commit-id, #mermaid-svg-KlQCVM1aHWdwKg9a .commit-msg, #mermaid-svg-KlQCVM1aHWdwKg9a .branch-label {fill: lightgrey;color: lightgrey;font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family); }#mermaid-svg-KlQCVM1aHWdwKg9a .pieTitleText {text-anchor: middle;font-size: 25px;fill: black;font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family); }#mermaid-svg-KlQCVM1aHWdwKg9a .slice {font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family); }#mermaid-svg-KlQCVM1aHWdwKg9a g.stateGroup text {fill: #9370DB;stroke: none;font-size: 10px;font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family); }#mermaid-svg-KlQCVM1aHWdwKg9a g.stateGroup text {fill: #9370DB;fill: #333;stroke: none;font-size: 10px; }#mermaid-svg-KlQCVM1aHWdwKg9a g.statediagram-cluster .cluster-label text {fill: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a g.stateGroup .state-title {font-weight: bolder;fill: black; }#mermaid-svg-KlQCVM1aHWdwKg9a g.stateGroup rect {fill: #ECECFF;stroke: #9370DB; }#mermaid-svg-KlQCVM1aHWdwKg9a g.stateGroup line {stroke: #9370DB;stroke-width: 1; }#mermaid-svg-KlQCVM1aHWdwKg9a .transition {stroke: #9370DB;stroke-width: 1;fill: none; }#mermaid-svg-KlQCVM1aHWdwKg9a .stateGroup .composit {fill: white;border-bottom: 1px; }#mermaid-svg-KlQCVM1aHWdwKg9a .stateGroup .alt-composit {fill: #e0e0e0;border-bottom: 1px; }#mermaid-svg-KlQCVM1aHWdwKg9a .state-note {stroke: #aaaa33;fill: #fff5ad; } #mermaid-svg-KlQCVM1aHWdwKg9a .state-note text {fill: black;stroke: none;font-size: 10px; }#mermaid-svg-KlQCVM1aHWdwKg9a .stateLabel .box {stroke: none;stroke-width: 0;fill: #ECECFF;opacity: 0.7; }#mermaid-svg-KlQCVM1aHWdwKg9a .edgeLabel text {fill: #333; }#mermaid-svg-KlQCVM1aHWdwKg9a .stateLabel text {fill: black;font-size: 10px;font-weight: bold;font-family: 'trebuchet ms', verdana, arial;font-family: var(--mermaid-font-family); }#mermaid-svg-KlQCVM1aHWdwKg9a .node circle.state-start {fill: black;stroke: black; }#mermaid-svg-KlQCVM1aHWdwKg9a .node circle.state-end {fill: black;stroke: white;stroke-width: 1.5; }#mermaid-svg-KlQCVM1aHWdwKg9a #statediagram-barbEnd {fill: #9370DB; }#mermaid-svg-KlQCVM1aHWdwKg9a .statediagram-cluster rect {fill: #ECECFF;stroke: #9370DB;stroke-width: 1px; }#mermaid-svg-KlQCVM1aHWdwKg9a .statediagram-cluster rect.outer {rx: 5px;ry: 5px; }#mermaid-svg-KlQCVM1aHWdwKg9a .statediagram-state .divider {stroke: #9370DB; }#mermaid-svg-KlQCVM1aHWdwKg9a .statediagram-state .title-state {rx: 5px;ry: 5px; }#mermaid-svg-KlQCVM1aHWdwKg9a .statediagram-cluster.statediagram-cluster .inner {fill: white; }#mermaid-svg-KlQCVM1aHWdwKg9a .statediagram-cluster.statediagram-cluster-alt .inner {fill: #e0e0e0; }#mermaid-svg-KlQCVM1aHWdwKg9a .statediagram-cluster .inner {rx: 0;ry: 0; }#mermaid-svg-KlQCVM1aHWdwKg9a .statediagram-state rect.basic {rx: 5px;ry: 5px; }#mermaid-svg-KlQCVM1aHWdwKg9a .statediagram-state rect.divider {stroke-dasharray: 10,10;fill: #efefef; }#mermaid-svg-KlQCVM1aHWdwKg9a .note-edge {stroke-dasharray: 5; }#mermaid-svg-KlQCVM1aHWdwKg9a .statediagram-note rect {fill: #fff5ad;stroke: #aaaa33;stroke-width: 1px;rx: 0;ry: 0; }:root {--mermaid-font-family: '"trebuchet ms", verdana, arial';--mermaid-font-family: "Comic Sans MS", "Comic Sans", cursive; }#mermaid-svg-KlQCVM1aHWdwKg9a .error-icon {fill: #552222; }#mermaid-svg-KlQCVM1aHWdwKg9a .error-text {fill: #552222;stroke: #552222; }#mermaid-svg-KlQCVM1aHWdwKg9a .edge-thickness-normal {stroke-width: 2px; }#mermaid-svg-KlQCVM1aHWdwKg9a .edge-thickness-thick {stroke-width: 3.5px; }#mermaid-svg-KlQCVM1aHWdwKg9a .edge-pattern-solid {stroke-dasharray: 0; }#mermaid-svg-KlQCVM1aHWdwKg9a .edge-pattern-dashed {stroke-dasharray: 3; }#mermaid-svg-KlQCVM1aHWdwKg9a .edge-pattern-dotted {stroke-dasharray: 2; }#mermaid-svg-KlQCVM1aHWdwKg9a .marker {fill: #333333; }#mermaid-svg-KlQCVM1aHWdwKg9a .marker.cross {stroke: #333333; }:root { --mermaid-font-family: "trebuchet ms", verdana, arial;}#mermaid-svg-KlQCVM1aHWdwKg9a {color: rgba(0, 0, 0, 0.75);font: ;}

生成随机系数矩阵A和常数项b
对数据进行截断使得两种程序中数据完全相等
分别用Matlab和Python求解
计算各自误差值
  • 关键问题一:保证两种语言处理的数据完全一致

这里第二步其实很重要,由于Matlab默认数据类型为double,而Python的默认类型为float,如果直接粘贴或者导入就有可能在Python中被截断,从而导致在一开始使用的数据就存在误差,使得实验结果不可信。

  • 关键问题二:选取合适的指标来进行对比

第四步主要是计算两种程序跑完后均方根误差:
∥Ax−b∥\|Ax-b\| ∥Ax−b∥
这种验证方式比较简单。之所以这么做是因为数值计算过程中不一定能够得到100%的精度,这个问题的主要原因是在数值计算时数字不会一直保持为真实值。尤其在除法时由于数值长度的限制会被强行截断一部分,因此用这个指标来验证方法是否正确是很有必要的。当然,根据它的定义不难看出,这个值最小精度就越高。

2、数据生成和处理

首先由Matlab随机生成一个矩阵和一个向量:

A = rand(5,5)
b = rand(5,1)

然后在Matlab中将数据粘至txt文件,再整理成Python的格式粘生成Numpy数组,同时也将数据重新粘回Matlab代码。这里原始数据是这个样子:

A的值

0.814723686393179    0.0975404049994095  0.157613081677548   0.141886338627215   0.655740699156587
0.905791937075619   0.278498218867048   0.970592781760616   0.421761282626275   0.0357116785741896
0.126986816293506   0.546881519204984   0.957166948242946   0.915735525189067   0.849129305868777
0.913375856139019   0.957506835434298   0.485375648722841   0.792207329559554   0.933993247757551
0.632359246225410   0.964888535199277   0.800280468888800   0.959492426392903   0.678735154857774

b的值

0.757700000000000
0.743100000000000
0.392200000000000
0.655500000000000
0.171200000000000

当然也可以采用round函数一类的方式进行处理。

注意 这种操作虽然看似简单并且有点笨拙,但能很好地保证二者在初始化时得到的结果是一样的。

Matlab的操作相对简单,直接将txt中的内容首尾加上[],再粘回代码即可。对于Python要稍加处理,最终的代码是:

import numpy as np
A = np.array([0.814723686393179,0.905791937075619,0.126986816293506,0.913375856139019,0.632359246225410,0.0975404049994095,0.278498218867048,0.546881519204984,0.957506835434298,0.964888535199277,0.157613081677548,0.970592781760616,0.957166948242946,0.485375648722841,0.800280468888800,0.141886338627215,0.421761282626275,0.915735525189067,0.792207329559554,0.959492426392903,0.655740699156587,0.0357116785741896,0.849129305868777,0.933993247757551,0.678735154857774]).reshape([5,5]).Tb = np.array([0.757700000000000,0.743100000000000,0.392200000000000,0.655500000000000,0.171200000000000])

注意这里Python的代码最后要进行一下转置,否则与Matlab的数据就就不一致。

整理完成看看数据和txt中的形状就完全一致了:

Aarray([[0.81472369, 0.0975404 , 0.15761308, 0.14188634, 0.6557407 ],[0.90579194, 0.27849822, 0.97059278, 0.42176128, 0.03571168],[0.12698682, 0.54688152, 0.95716695, 0.91573553, 0.84912931],[0.91337586, 0.95750684, 0.48537565, 0.79220733, 0.93399325],[0.63235925, 0.96488854, 0.80028047, 0.95949243, 0.67873515]])barray([0.7577, 0.7431, 0.3922, 0.6555, 0.1712])

这里b变成了一维数组但不影响计算。

3、求解方程组
  • Matlab部分

用Matlab求解线性方程组十分简单,直接一个反除号就可以搞定:

>> x = A \ bx =-0.84093.77013.8572-8.00502.4445

注意这里x跑出来的值只是显示值,实际值的位数是远大于这个值的:

-0.840928749143009   3.77006972003583    3.85716025378275    -8.00495386938441   2.44448015270832
  • Python部分

用Python求解线性方程组的方式其实也很多,最简单的办法就是直接用scipy库的linalg模块中的solve方法:

from scipy.linalg import solvex = solve(A,b)xarray([-0.84092875,  3.77006972,  3.85716025, -8.00495387,  2.44448015])

为了更准确地对比数值的精度,我们将x的值用np.savetxt函数保存至文件再看具体数字:

-8.409287491430091910e-01
3.770069720035832184e+00
3.857160253782746739e+00
-8.004953869384410226e+00
2.444480152708315313e+00

仔细观察发现,Python输出的数值位数要稍多一些,这是因为我们用Matlab粘贴的时候它默认显示的位数不太一样。

注:这也是为什么我们一定要在最开始进行数据统一规范的原因。Matlab在存储变量时数据是保存了完整的位数的,但即便从变量查看窗口粘出来的数值也不一定是完整的数字。

4、精度验证
  • Matlab部分
error = norm(A*x-b)error =1.5801e-15% 查看15位
% 1.580118849929844e-15
  • Python部分
np.sqrt(np.sum((A.dot(x)-b)**2))1.580118849929844e-15

注意,由于两个数的量级都是1e-15,因此这里只需要看它们各自前16位即可。仔细观察不难发现,两个程序得到的结果完全一致。

这里我们把有效数字从第一位到最后一位排列直接放出来对比

# python
1580118849929843973242195223453751011511480618378089335607228349545039236545562744140625
% matlab
1580118849929843973242195223453751011511480618378089335607228349545039236545562744140625

肉眼观察信不过,全部放入Python字符串试试对不对:

#第一行数字(python结果):
s1 = '1580118849929843973242195223453751011511480618378089335607228349545039236545562744140625'#第二行数字(matlab结果):
s2 = '1580118849929843973242195223453751011511480618378089335607228349545039236545562744140625'

判断一下是否相等:

s1 == s2True

搞定!

5、结论

Python(Numpy, Scipy)的数值精度与Matlab是完全一致的。

其本质的原因来自于二者在矩阵计算都是使用的LAPACK, 详情可见:
Matlab和Python(Numpy,Scipy)与Lapack的关系

深度对比Python(Numpy,Scipy)与Matlab的数值精度相关推荐

  1. matlab error函数_深度对比Python(Numpy,Scipy)与Matlab的数值精度

    CSDN原文:https://blog.csdn.net/cauchy7203/article/details/107785295 Matlab一度被认为是最专业的数值计算工具之一,相信许多同学都或多 ...

  2. Matlab和Python(Numpy,Scipy)与Lapack的关系

    说到数值计算,可能许多人都能立马想到Matlab.Matlab多年的持续影响力已经让它成为许多人心中科学计算的代名词.但它底层一个重要的库Lapack却很少有人知道. 而Python年龄比Matlab ...

  3. python用于统计学_R 和 Python (numpy scipy pandas) 用于统计学分析,哪个更好?

    用R做过Python擅长的,也用Python做过R擅长的.有一天,打开官网,看到两句话,豁然开朗. R:R is a free software environment for statistical ...

  4. python apriori算法 sklearn_R 和 Python (numpy scipy pandas) 用于统计学分析,哪个更好?...

    可以两个一起学,参见我的博文,可以做个索引. Python 和 R 数据分析/挖掘工具互查 写在前面 在此总结一些在数据分析/挖掘中可能用到的功能,方便大家索引或者从一种语言迁移到另一种.当然,这篇博 ...

  5. 使用python+numpy+scipy进行图像处理实战

    以前照相没有像现在这样那么容易的,而在现在你只需要一部手机,就可以免费拍照,而在上一代人之前,业余艺术家和真正的艺术家拍照的费用非常昂贵,并且每张照片的成本也不是免费的. 我们拍照是为了及时地保存美好 ...

  6. 深度对比Python的3种“字符串格式化”方法,看看你喜欢哪一种?

    3种字符串格式化工具的简单介绍 python2.5版本之前,我们使用的是老式字符串格式化输出%s. 从python3.0版本开始起(python2.6同期发布),Python中同时支持两个版本的格式化 ...

  7. 深度对比 Python 与 Java 的区别(一)

    引入 高中有一好友,在大学期间苦练 Java,各类八股文烂熟于心,最终进入某大厂却在维护 Python 项目. 而本人不思进取,不想背八股文,于是大学期间只是苟着写 Python,然而却最终进入某互联 ...

  8. python基础知识及数据分析工具安装及简单使用(Numpy/Scipy/Matplotlib/Pandas/StatsModels/Scikit-Learn/Keras/Gensim))

    Python介绍. Unix & Linux & Window & Mac 平台安装更新 Python3 及VSCode下Python环境配置配置 python基础知识及数据分 ...

  9. python的数据与matlab互通:SciPy

    有时候需要利用python进行科学计算,但需要Matlab进行交互式画图,因此需要掌握python与matlab数据互通的技巧:SciPy SciPy 提供了与 Matlab 的交互的方法. SciP ...

最新文章

  1. steam自建服务器游戏_虽有差评销量却还是直步青云,《Atlas》力登Steam榜单前茅...
  2. 微信小程序打开PDF
  3. SharePoint designer 文件--新建中没有工作流
  4. 从前有座山,山里有座庙:递归之法
  5. 如何对数据库中的表以及表中的字段进行重命名
  6. [css] box-sizing的宽度包含了哪些?
  7. git 获取最新代码_程序员必知:这是一份全面 amp; 详细的 Git与Github 介绍指南
  8. 7. vue-cli 安装和使用脚手架
  9. 自定义view圆环的改变
  10. jQuery框架学习
  11. C语言解释器的实现--让脚本跑起来(六)
  12. 流光快门Matlab,手机相机中的流光快门怎么用?教你用流光快门拍出最炫酷的照片...
  13. SecureCRT and SecureFX 8.3 中文版
  14. 好用过头的LeetCode刷题模板分享!(已拿亚麻offer)
  15. matlab灰色预测关联度
  16. [VS2010]逸雨清风 永久稳定音乐外链生成软件V0.1
  17. 《测绘管理与法律法规》真题易错本
  18. NS3仿真之LTE数据分析RSRP,SINR,吞吐量
  19. ES学习看这一篇文章就够了
  20. Proteus仿真最小系统板的绘制及流水灯

热门文章

  1. 检查卷位图时发现损坏怎么修复_中频弯管严密性如何测试?怎么修复中频弯管?...
  2. linux进程假死的原因_一次Spring Boot假死诊断
  3. Windows server2008修改远程桌面端口的方法
  4. 问题集合---《平时遇到的问题 + 参考解决方式》
  5. 解决AttributeError: 'module' object has no attribute 'main' 安装第三方包报错
  6. Document 对象描述
  7. css实现倒8字效果
  8. 关于ios app发布的中间证书的要求--解决WWDR证书过期方案
  9. 摩托罗拉SE4500 德州仪器TI Omap37xx/AM3715/DM3730/AM3530 wince6.0/Windows Mobile 6.5平台 二维软解调试记录及相关解释
  10. WINCE cvrtbin命令简介