CN114091363A - 基于量子算法的计算流体动力学模拟方法、装置及设备 - Google Patents

基于量子算法的计算流体动力学模拟方法、装置及设备 Download PDF

Info

Publication number
CN114091363A
CN114091363A CN202010771592.2A CN202010771592A CN114091363A CN 114091363 A CN114091363 A CN 114091363A CN 202010771592 A CN202010771592 A CN 202010771592A CN 114091363 A CN114091363 A CN 114091363A
Authority
CN
China
Prior art keywords
quantum
oracle
state
grid
line
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202010771592.2A
Other languages
English (en)
Other versions
CN114091363B (zh
Inventor
李叶
安宁波
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Origin Quantum Computing Technology Co Ltd
Original Assignee
Origin Quantum Computing Technology Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Origin Quantum Computing Technology Co Ltd filed Critical Origin Quantum Computing Technology Co Ltd
Priority to CN202010771592.2A priority Critical patent/CN114091363B/zh
Priority to US18/017,571 priority patent/US20230297745A1/en
Priority to PCT/CN2020/140833 priority patent/WO2022027916A1/zh
Priority to EP20948742.0A priority patent/EP4195090A4/en
Publication of CN114091363A publication Critical patent/CN114091363A/zh
Application granted granted Critical
Publication of CN114091363B publication Critical patent/CN114091363B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06NCOMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N10/00Quantum computing, i.e. information processing based on quantum-mechanical phenomena
    • G06N10/20Models of quantum computing, e.g. quantum circuits or universal quantum computers
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • G06F30/28Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/10Geometric CAD
    • G06F30/15Vehicle, aircraft or watercraft design
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • G06F30/23Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06NCOMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N10/00Quantum computing, i.e. information processing based on quantum-mechanical phenomena
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06NCOMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N10/00Quantum computing, i.e. information processing based on quantum-mechanical phenomena
    • G06N10/60Quantum algorithms, e.g. based on quantum optimisation, quantum Fourier or Hadamard transforms
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00Details relating to CAD techniques
    • G06F2111/10Numerical modelling
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00Details relating to the application field
    • G06F2113/08Fluids
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14Force analysis or force optimisation, e.g. static or dynamic forces
    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02TCLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • Evolutionary Computation (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Computer Hardware Design (AREA)
  • Computing Systems (AREA)
  • Mathematical Physics (AREA)
  • Computational Mathematics (AREA)
  • Algebra (AREA)
  • Fluid Mechanics (AREA)
  • Artificial Intelligence (AREA)
  • Condensed Matter Physics & Semiconductors (AREA)
  • Data Mining & Analysis (AREA)
  • Software Systems (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Automation & Control Theory (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本申请实施例提供了基于量子算法的计算流体动力学模拟方法、装置及设备,量子计算的固有特征是量子并行性。在希尔伯特空间中,量子比特的量子态是一个向量,每个作用在量子比特上的量子逻辑门操作都是在希尔伯特空间中的量子态变换。如果将对一个实际问题的描述编码为一组量子比特的量子状态向量,则可以通过执行操作来改变向量实现问题的求解,而无需迭代每个量子比特的状态。例如,如果将矩阵和向量相乘,不再需要for循环,仅需要量子逻辑门操作作用于输入的量子态上即可。本申请实施例中,实现了利用量子计算实现CFD模拟,利用量子计算,相对于经典算法有指数级的加速,从而能够降低CFD的复杂度,增加CFD的实用性。

Description

基于量子算法的计算流体动力学模拟方法、装置及设备
技术领域
本申请涉及量子计算技术领域,特别是涉及基于量子算法的计算流体动力学模拟方法、装置及设备。
背景技术
CFD(Computational Fluid Dynamics,计算流体动力学)是近代流体力学,数值数学和计算机科学结合的产物,是一门具有强大生命力的交叉科学。它以电子计算机为工具,应用各种离散化的数学方法,对流体力学的各类问题进行数值实验、计算机模拟和分析研究,以解决各种实际问题。CFD是模拟飞机气动设计的一个基本工具,从中可以获得飞机周围的流动特性和作用于飞机表面的气动力。这些计算在计算上是昂贵的,这使得寻求高效的CFD算法成为研究人员永远没有解决的问题。
现有技术中,在通过计算流体动力学模拟流体运动的过程中,针对隐式方案是无条件稳定的,需要求解线性系统
Figure BDA0002616833180000011
其中
Figure BDA0002616833180000012
是一个系数矩阵,在计算过程中必须求逆。对该矩阵求逆的计算非常复杂,并且需要采用迭代方法。该系数矩阵的维数跟网格节点成线性关系,当网格节点规模很大时,使用例如,Jacobi(雅可比)方法,LU方法或GMRES(广义最小残差法)等经典的算法,实现该系数矩阵的求解将非常的困难。
量子计算的固有特征是量子并行性。在希尔伯特空间(复线性空间)中,N量子比特的量子态是一个向量,每个量子操作都是在此线性空间中的变换。如果将对一个实际问题的描述编码为一组量子比特的量子状态向量,则可以通过并行执行一个操作来处理向量实现问题的求解,而无需迭代一组量子比特中每个量子比特的状态。例如,如果将矩阵和向量相乘,不再需要for循环,仅需要几个量子门作用于输入的量子态上即可。因此量子计算的并行性可以大大增加CFD的计算速度。所以希望利用量子计算实现CFD模拟。
发明内容
本申请实施例的目的在于提供一种基于量子算法的计算流体动力学模拟方法、装置及设备,以实现利用量子计算实现CFD模拟。具体技术方案如下:
第一方面,本申请实施例提供了一种基于量子算法的计算流体动力学模拟方法,所述方法包括:
在利用有限体积法解析计算流体动力学的过程中,针对流体运动的离散化数值网格中各网格单元,构建表示所述网格单元坐标信息的第一量子线路、表示所述网格单元的状态参数的第二量子线路,其中,所述网格单元的状态参数存储在量子随机存取存储器中,所述量子随机存取存储器可操作处于量子叠加态的地址和数据;
基于所述第一量子线路、所述第二量子线路、所述量子随机存取存储器构建表示所述网格单元流体状态变化的线性系统方程的参数的第三量子线路;
针对所有的网格单元,基于所述第三量子线路求解所述网格单元的线性系统方程,当所述网格单元的流体状态趋于稳定时,获得所述网格单元线性系统方程的状态参数表示的流体状态作为所述网格单元的目标状态。
第二方面,本申请实施例提供了一种基于量子算法的计算流体动力学模拟装置,所述装置包括:
网格量子线路构建模块,用于在利用有限体积法解析计算流体动力学的过程中,针对流体运动的离散化数值网格中各网格单元,构建表示所述网格单元坐标信息的第一量子线路、表示所述网格单元的状态参数的第二量子线路,其中,所述网格单元的状态参数存储在量子随机存取存储器中,所述量子随机存取存储器可操作处于量子叠加态的地址和数据;
参数量子线路构建模块,用于基于所述第一量子线路、所述第二量子线路、所述量子随机存取存储器构建表示所述网格单元流体状态变化的线性系统方程的参数的第三量子线路;
目标状态获取模块,用于针对所有的网格单元,基于所述第三量子线路求解所述网格单元的线性系统方程,当所述网格单元的流体状态趋于稳定时,获得所述网格单元线性系统方程的状态参数表示的流体状态作为所述网格单元的目标状态。
第三方面,本申请实施例提供了一种量子计算机设备,包括:量子线路、量子效应器及量子随机存取存储器,所述量子计算机设备在运行时实现任一所述的基于量子算法的计算流体动力学模拟方法。
第四方面,本申请实施例提供了一种计算机可读存储介质,其特征在于,所述计算机可读存储介质内存储有计算机程序,所述计算机程序被处理器执行时实现任一所述的基于量子算法的计算流体动力学模拟方法。
本申请实施例提供的基于量子算法的计算流体动力学模拟方法、装置及设备,在利用有限体积法解析计算流体动力学的过程中,针对流体运动的离散化数值网格中各网格单元,构建表示网格单元坐标信息的第一量子线路、表示网格单元的状态参数的第二量子线路,其中,网格单元的状态参数存储在量子随机存取存储器中,量子随机存取存储器可操作处于量子叠加态的地址和数据;基于第一量子线路、第二量子线路、量子随机存取存储器构建表示网格单元流体状态变化的线性系统方程的参数的第三量子线路;针对所有的网格单元,基于第三量子线路求解网格单元的线性系统方程,当网格单元的流体状态趋于稳定时,获得网格单元线性系统方程的状态参数表示的流体状态作为网格单元的目标状态,实现了利用量子计算实现CFD模拟。量子计算的固有特征是量子并行性,利用量子计算,相对于经典算法有指数级的加速,从而能够降低CFD的复杂度,增加CFD的实用性。当然,实施本申请的任一产品或方法并不一定需要同时达到以上所述的所有优点。
附图说明
为了更清楚地说明本申请实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本申请的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1为相关技术中离散化数值网格的一种示意图;
图2为相关技术中离散化数值网格的另一种示意图;
图3为相关技术中迭代时间与网格单元个数关系的曲线图;
图4本申请实施例的基于量子算法的计算流体动力学模拟方法的第一种示意图;
图5为本申请实施例的初始状态的量子线路的一种示意图;
图6为本申请实施例的数字版本的
Figure BDA0002616833180000041
的量子线路的一种示意图;
图7为本申请实施例的Oracle Ob的量子线路的一种示意图;
图8为本申请实施例的Oracle Ol的量子线路的一种示意图;
图9为本申请实施例的Oracle OA的量子线路的一种示意图;
图10为本申请实施例的量子线性求解器的量子线路的一种示意图;
图11为本申请实施例的后处理过程中的量子线路的一种示意图;
图12a为现有技术中经典模型的一种示意图;
图12b为本申请实施例的量子-经典混合模型的一种示意图;
图13为本申请实施例的基于量子算法的计算流体动力学模拟方法的第二种示意图;
图14为本申请实施例的网格单元个数与量子比特数之间关系的曲线图;
图15为本申请实施例的QDAC算法的量子线路的一种示意图;
图16a为本申请实施例的实部QADC算法的量子线路的一种示意图;
图16b为本申请实施例的相位估计的量子线路的一种示意图;
图16c为本申请实施例的实部QADC的量子线路Q部分的量子线路的一种示意图;
图17a为本申请实施例的关于T的量子线路的一种示意图;
图17b为本申请实施例的关于K的量子线路的一种示意图;
图17c为本申请实施例的关于W的量子线路的一种示意图;
图18为本申请实施例的基于量子算法的计算流体动力学模拟装置的一种示意图。
具体实施方式
下面将结合本申请实施例中的附图,对本申请实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本申请一部分实施例,而不是全部的实施例。基于本申请中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本申请保护的范围。
申请概述:
CFD模拟是飞机气动设计的一个基本工具,从中可以获得飞机周围的流动特性和作用于飞机表面的气动力。CFD问题通常与用数值方法求解一系列方程有关。其中,CFD求解器需要对偏导数方程进行时间积分,因此需要使用时间增量进行迭代以计算系统的小变化,并在系统变得稳定时停止。在这个问题设置上,经典算法中一般通过SU2(TheStandford University Unstructured suite,一种斯坦福大学开发的高精度偏微分方程求解器)进行计算。在CFD模拟时,从数学模型开始,需要定义以下元素,包括:数值网格、离散化方法、数值方案及解决方法。
数值网格:网格是解决问题的空间连续域的离散表示。将每个空间点都认为是一个节点,由一组相邻节点界定的体积定义一个网格单元。在结构化网格的情况下,每个节点和网格单元由一组索引唯一标识,例如,在二维情况下可以表示为(a,b),其中,a与b分别表示二维空间的横坐标及纵坐标,例如图1所示,由一组网格节点(A,B,C,D,…)组成已知坐标的网格,控制体积中心(又叫单元)坐标为P(a,b)。非结构化网格通常由四面体元素组成,并且需要一个连通性列表,该列表指定给定的一组节点组成一个单元的方式。
离散化方法:解决上述流体动力学方程式的标准方法是FVM(Finite volumemethod,有限体积法)。FVM方法中将数值网格划分为多个网格单元,考虑了与每个网格单元相关联的局部体积,并应用了积分守恒定律。这种定律表明,体积Ω内流体的所有变量U的变化仅取决于其表面S上的通量F。FVM可以自动满足流体力学方程式的守恒性质。这个守恒方程可以写成如下:
Figure BDA0002616833180000051
数值方案:在空间连续域离散化之后,要求解的方程的数学运算符也必须离散化,方程的这种近似是通过数值方案实现的。考虑两种类型的数值方案:时间积分方案和空间方案。时间积分方案是指从流体的最初状态开始对需要求解的变量的时间进行积分建模,而空间方案则代表这些需要求解的变量的空间梯度。
解决方法:一旦建立了离散问题和数值方案,就必须解决一个代数方程组。根据方案的性质,需要迭代方法。
求解公式(1)守恒方程的一种有效方法是FVM。基于上述介绍的离散化数值方案和解决方法,需要对空间的数值网格进行求解。在二维空间,每个空间点是一个节点,由一组索引标识。通过最近相邻节点连接的节点集合形成网格,示例性的例如图2所示,节点(a,b)、节点(a,b+1)、节点(a+1,b)、节点(a+1,b+1)组成的四边形网格。两个相邻节点之间的线的中点和相邻网格的中心点决定了一个节点周围的单元,例如图2中的以节点(a,b)为中心的八边形单元。守恒方程应用于网格的每个单元,即网格单元。在公式(1)的离散化表示中,微分运算被数值逼近代替,并选择隐式欧拉格式来处理时间微分。公式(1)离散为:
Figure BDA0002616833180000064
Un为第n个时刻(第n次积分)所有网格单元的流体状态,Un={Ui},i=0,1,2,…,N-1,Ui=(ρ,ρu,ρν,ρE)i表示第i个网格单元的流体状态,ρ表示流体密度,u及v表示流体速度在二维坐标下的两个分量,E表示单位质量流体的总能量,bn(Un)表示第n个时刻的残差量,且其为与Un相关的函数。
令ΔUn=Un+1-Un,则:
Figure BDA0002616833180000061
将公式(3)代入公式(2),得到:
Figure BDA0002616833180000062
其中,δi,i′δj,j′为与网格单元的坐标信息相关的参数。
Figure BDA0002616833180000063
则,公式(4)等价为:
A(Un)ΔUn=-bn(Un) (5)
其中,A(Un)表示第n个时刻的系数矩阵,且其为与Un相关的函数。
将A(Un)简写为An,将-bn(Un)简写为bn,从而得到了一个需要求解线性系统:AnΔUn=bn的迭代问题。
当网格单元数目为N时,那么Un的维数为4N,求解时,获取这N个网络单元的流体初始状态得到U0,根据U0计算得到A0,b0,解A0(U1-U0)=b0,得到U1,然后计算A1,b1,得到U2,一直迭代下去,迭代截止条件是bn趋近于0。在这种情况下,需要使用求解AnΔUn=bn形式的线性系统的方法,其中,ΔUn=Un+1-Un
经典方案中,隐式欧拉求解器具有线性复杂度,An是一个稀疏矩阵,在计算过程中必须求逆。对该矩阵求逆的计算非常复杂,并且需要采用迭代方法。经典的求解方法采用多重网格法,该稀疏矩阵的维数跟网格节点成线性关系,复杂度为O(N),其中N为网格单元数。经典方案中,运行时间与网格单元的数量具有很好的线性关系,例如图3,当网格节点规模很大时,使用例如,Jacobi(雅可比)方法,LU方法或GMRES(广义最小残差法)等经典的算法,实现将非常的困难。
有鉴于此,为了降低CFD的复杂度,增加CFD的实用性,本申请实施例提供了一种基于量子算法的计算流体动力学模拟方法,该方法包括:在利用有限体积法解析计算流体动力学的过程中,针对流体运动的离散化数值网格中各网格单元,构建表示所述网格单元坐标信息的第一量子线路、表示所述网格单元的状态参数的第二量子线路,其中,所述网格单元的状态参数存储在量子随机存取存储器中,所述量子随机存取存储器可操作处于量子叠加态的地址和数据;基于所述第一量子线路、所述第二量子线路、所述量子随机存取存储器构建表示所述网格单元流体状态变化的线性系统方程的参数的第三量子线路;针对所有的网格单元,基于所述第三量子线路求解所述网格单元的线性系统方程,当所述网格单元的流体状态趋于稳定时,获得所述网格单元线性系统方程的状态参数表示的流体状态作为所述网格单元的目标状态。
量子计算的固有特征是量子并行性。在希尔伯特空间(复线性空间)中,N量子比特的量子态是一个向量,每个量子操作都是在此线性空间中的变换。如果将对一个实际问题的描述编码为一组量子比特的量子状态向量,则可以通过并行执行一个操作来处理向量实现问题的求解,而无需迭代一组量子比特中每个量子比特的状态。例如,如果将矩阵和向量相乘,不再需要for循环,仅需要几个量子门作用于输入的量子态上即可。本申请实施例中,实现了利用量子计算实现CFD模拟,利用量子计算,相对于经典算法有指数级的加速,从而能够降低CFD的复杂度,减少CFD的计算时间,增加CFD的实用性。
以下进行具体说明:
首先,对本申请实施例中的符号进行解释。
N:网格单元数。
Figure BDA0002616833180000081
表示编码所有网格单元初始状态所用的量子比特数。
qf:用浮点数表示的量子比特数。
ITimes:迭代次数。
Ui=(ρ,ρu,ρν,ρE)i,第i个网格单元的变量(也称为状态参数),也表示i个网格单元的流体状态,共有4个变量。
Un={Ui},i=0,1,2,…,N-1,Un表示问题对应的所有变量,也表示所有网格单元的流体状态的集合,维数是4N。
当线性系统方程AnΔUn=bn的残差bn趋于0时,可以得到稳定的解,量子算法是求解此类线性系统的有效方法。可以建立量子线路来完成计算。参见图4,本申请实施例提供了一种基于量子算法的计算流体动力学模拟方法,包括:
S11,在利用有限体积法解析计算流体动力学的过程中,针对流体运动的离散化数值网格中各网格单元,构建表示上述网格单元坐标信息的第一量子线路、表示上述网格单元的状态参数的第二量子线路,其中,上述网格单元的状态参数存储在量子随机存取存储器中,上述量子随机存取存储器可操作处于量子叠加态的地址和数据。
本申请实施例的基于量子算法的计算流体动力学模拟方法,可以通过具备量子计算能力的电子设备,例如量子计算机等实现,也可通过经典计算机模拟量子计算实现。
第一量子线路用于表示网格单元的坐标信息,第二量子线路用于表示网格单元的状态参数。例如针对第i个网格单元(或称为网格单元i):Ui=(ρ,ρu,ρν,ρE)i,第一量子线路用于表示i相关的坐标信息,第二量子线路用于表示(ρ,ρu,ρν,ρE)。
网格单元的状态参数存储在QRAM(Quantum Random Access Memory,量子随机存取存储器中,QRAM是经典随机存取存储器的量子模拟,是一种可以有效地存储和加载量子信息的结构。与普通方式不同,量子随机存取存储器可操作处于量子叠加态的地址和数据,QRAM在叠加中操作地址。当地址处于叠加状态时,QRAM使用户能够有效地获取QRAM中的数据。在一种可能的实施方式中,上述量子随机存取存储器内存储有的地址和数据的预设映射关系,上述地址用第一叠加态编码,上述第一叠加态中的每个本征态对应一个地址信息,上述数据用第二叠加态编码,上述第二叠加态中的每个本征态对应一个数据信息。本申请实施例中QRAM既可以采用量子的方式实现,也可以采用经典计算机模拟实现,利用经典计算机模拟实现QRAM为现有技术,还可以采用硬件电路实现,此处不再赘述,
S12,基于上述第一量子线路、上述第二量子线路、上述量子随机存取存储器构建表示上述网格单元流体状态变化的线性系统方程的参数的第三量子线路。
利用第一量子线路可以得到各网格单元的坐标信息,利用第二量子线路可以得到各网格单元的初始状态参数,由因为网格单元的流体状态时通过状态参数表示的,因此可以得到各网格单元的初始流体状态。量子随机存取存储器能够存储各网格单元的状态参数(包括初始状态参数及迭代后的状态参数),因此基于第一量子线路、第二量子线路、量子随机存取存储器,可以构建表示网格单元流体状态变化的线性系统方程的参数的第三量子线路。
S13,针对所有的网格单元,基于上述第三量子线路求解上述网格单元的线性系统方程,当上述网格单元的流体状态趋于稳定时,获得上述网格单元线性系统方程的状态参数表示的流体状态作为上述网格单元的目标状态。
可以通过第三量子线路对所有的网格单元的流体状态进行迭代,当网格单元的流体状态趋于稳定时,即残差bn趋于0时,可以得到此时各网格单元线性系统方程的流体状态,即为各网格单元的目标状态。
可选的,在得到各网格单元的目标状态后,可以将网格单元的目标状态转化为模拟编码态,并测量上述模拟编码态,从而完成流体运动的计算流体动力学的解析。
从量子随机存取存储器中读取的网格单元的目标状态为数字态(量子数字态),需要转化为模拟编码态(量子模拟态),并通过测量模拟编码态得到用于描述流体运动的经典数据,从而完成后续流体运动的计算流体动力学的解析。
网格单元的目标状态的量子形式的数字编码态表示为
Figure BDA0002616833180000101
保存在QRAM中。经典信息的提取过程如下:首先从QRAM中读取数字编码态的
Figure BDA0002616833180000102
然后通过QDAC(QuantumDigital-to-Analog Converter,量子数模转换器),将数字编码态转化为模拟编码态:
Figure BDA0002616833180000103
其中,
Figure BDA0002616833180000104
可以理解的是,模拟编码态是指数字信息编码在量子态振幅上的态,而数字编码态是指数字信息编码在量子态本身(相等于寄存器)上的态。最后,通过测量模拟编码态
Figure BDA0002616833180000105
来提取经典信息。
在计算的过程中需要用到N个网格单元,当网格单元的数量N很大时,可选的,可以分多次采集网格单元的流体状态。单次采样时,可以只采样N个网格单元中的N1个网格单元的流体状态,从而通过多次采样得到最后输出的经典信息。假设降维后第j个网格包含aj个原来的网格,一共采样M次,第j个网格出现的次数为bj次,则得到的N1个网格的解可以表示为:
Figure BDA0002616833180000106
Figure BDA0002616833180000107
最终采样M次,由于QDAC是概率性的操作,因此有概率
Figure BDA0002616833180000108
会将数字信息成功转化为模拟信息,因此需要M/p次从QRAM中提取信息以及M/p次的QDAC转换。
量子计算的固有特征是量子并行性,本申请实施例中,利用量子计算的量子并行性,相对于经典算法有指数级的加速,从而能够降低CFD计算的复杂度,减少CFD计算的时间,增加CFD的实用性。
为了便于量子线路的设计,方便明确各量子线路的功能,在一种可能的实施方式中,上述第一量子线路包括第一Oracle及第二Oracle;上述第一Oracle用于提取指定网格单元的相邻网格数,上述第二Oracle用于提取网格单元的坐标信息。
利用网格化方案将网格单元坐标相关的参数编码到第一量子线路中。此处可以定义了两个Oracle:第一Oracle(以下表示为Oadjacent)和第二Oracle(以下表示为Ocoor)。其中,Oadjacent用于提取指定网格单元的相邻网格数,可以表示为:
Oadjacent|i>|j>=|i>|g(i,j)> (6)
其中:g(i,j)是网格单元i的第j个相邻的网格单元。
Ocoor用于提取网格单元的坐标,可以表示为:
Figure BDA0002616833180000111
其中
Figure BDA0002616833180000112
表示和网格单元i相关的坐标信息,比如网格单元i各个顶点的坐标。还可以用Ocoor得到其它的位置信息,比如网格单元i中心的坐标和网格单元i每条边的法线。
构建Oadjacent和Ocoor所需的量子比特数量分别为2(m-2),m-2+C0qf,这里C0qf是一个量子比特序列,表示从网格单元坐标中获得的坐标信息和其他位置信息,C0是一个常数。构建Oadjacent和Ocoor的复杂度依赖于网格方案的复杂度。其中,复杂度可以为
Figure BDA0002616833180000114
所以Oadjacent和Ocoor的复杂度为
Figure BDA0002616833180000113
在本申请实施例中,第一量子线路包括第一Oracle及第二Oracle;第一Oracle用于提取指定网格单元的相邻网格数,第二Oracle用于提取网格单元的坐标信息,将第一量子线路的功能分配给第一Oracle及第二Oracle实现,明确了每部分的功能,便于量子线路的设计及应用。
在数据转换的过程中,需要用到量子数字态及模拟编码态两种形式的状态参数,在一种可能的实施方式中,上述第二量子线路,包括:通过量子态的振幅提取上述网格单元的状态参数的第二模拟态量子线路,和通过量子态的本征态提取上述网格单元的状态参数的第二数字态量子线路。
第二量子线路包括第二模拟态量子线路及第二数字态量子线路,第二模拟态量子线路用于通过量子态的振幅提取上述网格单元的状态参数,即用于表示模拟编码态的状态参数;第二数字态量子线路用于通过量子态的本征态提取上述网格单元的状态参数,即用于表示数字态的状态参数。
在利用FVM求解的过程中,每个网格单元的状态参数Ui的初值可以认为是相同的,可以表示为Ui=(x0,x1,x2,x3),其中,x0,x1,x2,x3分别为ρ,ρu,ρν,ρE的具体数值,所以所有网格单元的初始状态U0可以表示为:
Figure BDA0002616833180000121
其中:j=imod4,即i/4取余数,最终xj是x0、x1、x2和x3中的一个。
Figure BDA0002616833180000122
为U0的模拟编码态的表示形式,为可以用如图5所示的量子线路来制备。
同样,还可以构建数字版本的
Figure BDA0002616833180000123
即:
Figure BDA0002616833180000124
其中,
Figure BDA0002616833180000125
为U0的数字态的表示形式,对应的量子线路可以如图6所示。
Figure BDA0002616833180000126
可以保存在QRAM,QRAM的定义如下:
Figure BDA0002616833180000127
上式中的MU作为保存网格单元参数状态变量的数据存储器,是一个量子随机存取存储器QRAM,具有量子逻辑门的功能,并且MU可以作用在一个叠加态上面。需要说明的是,上述定义的QRAM是一种可以有效地存储和加载量子信息的结构,QRAM将地址和数据均编码在量子比特上,可以通过量子比特的叠加态表示地址,因此能够在叠加态中操作地址。QRAM能够使数字形式的量子数据并行加载和存储以时间复杂度O(logN)进行。其中,时间复杂度O(logN)表示当数据增大N倍时,耗时增大logN倍,其中的log可以是以2为底的。
由上述描述可知,线性系统方程可以为:AnΔUn=bn,其参数为系数矩阵An和当前时刻的残差量bn,因此可以将第三量子线路的功能进行更加细致的划分。
在一种可能的实施方式中,上述线性系统方程为:AnΔUn=bn;其中:ΔUn=Un+1-Un;其中:n表示记录流体状态的第n个时刻;Un表示第n个时刻网格单元的流体状态,Un={Ui},i=0,1,2,…,N-1,Ui表示网格单元i的流体状态;An表示第n个时刻的系数矩阵,且
Figure BDA0002616833180000131
表示An为与
Figure BDA0002616833180000132
相关的函数,bn表示第n个时刻的残差量,且
Figure BDA0002616833180000133
表示bn为与
Figure BDA0002616833180000134
相关的函数,N为网格单元的总数量,
Figure BDA0002616833180000135
为上述网格单元的坐标信息。
上述线性系统方程的参数包括当前时刻的系数矩阵An和当前时刻的残差量bn;上述第三量子线路包括对应上述系数矩阵An的Oracle OA,和对应上述残差量bn的Oracle Ob;系数矩阵An是一个稀疏矩阵,并且当i和j相邻时,元素Ai′j′是非零的;上述Oracle OA用于提取当前n时刻系数矩阵An的元素Ai′j′,其中,Ai′j′表示系数矩阵中第i′行第个j′元素;An={Ai′j′},i′=4i+i1;j′=4j+j1,i,j=0,1,2,…,N-1;i1,j1∈[0,3]的整数;对网格单元i的第j个相邻网格单元,上述Oracle OA的作用为OA|i>|j>=|i>|j>|Ai′j′>;上述Oracle Ob用于提取当前n时刻的残差量bn的元素bi,其中,bi表示网格单元i的残差量,bn={bi},i=0,1,2,…,N-1;对网格单元i,上述Oracle Ob的作用为Ob|0>=∑ibi|i>。
Oracle Ob对网格单元的初始流体状态,即|0>态的作用为Ob|0>=∑ibi|i>,bi表示网格单元i的残差量,为第
Figure BDA0002616833180000136
个网格单元及其邻近网格单元的状态参数和坐标的函数,这些函数的表达式是规则的,具体来说,bi的表达式只与i mod4相关,因此可以有效地构造bi,其中,mod表示求余运算。Oracle Ob的一种可能的量子线路图可以如图7所示。
首先,得到相关的状态参数和坐标:|i>|Urelated>|Xrelated>,Xrelated表示相关网格单元的坐标信息,Urelated表示相关网格单元对应的状态参数信息。
然后,计算bi得到:|i>|Urelated>|Xrelated>|bi>。
此处,可以不计算|Urelated>和|Xrelated>,|i>|bi>,因为其是为了计算目标状态而构建出的辅助态,可以使用QDAC进行操作,完成∑ibi|i>的过程。
在这个过程中,所需的量子比特的数量是m+Cbqf,m为一个常数,表示备用的量子比特的数量,Cbqf是一个辅助量子比特序列用来计算bi;也可以用于QDAC量子线路。假设QDAC的成功率是pQDAC,则MU,Oadjacent和Ocoor的使用次数为O(1/pQDAC)。量子逻辑门的复杂度是O(poly(m)/pQDAC+poly(1/∈pQDAC)),其中∈代表QDAC的准确度。
下面对Oracle Ob的一种可能的构建方式进行说明,在一种可能的实施方式中,上述基于上述第一量子线路、上述第二量子线路、上述量子随机存取存储器构建表示上述网格单元流体状态变化的线性系统方程的参数的第三量子线路,包括:
步骤A,针对bn包括的对应于网格单元i的残差量bi,获取网格单元i的相关网格单元。
步骤B,基于上述第二Oracle获取上述相关网格单元的坐标信息作为第一相关坐标信息Xrelated
步骤C,基于上述第二量子线路和上述量子随机存取存储器,获取相关网格单元对应的状态参数信息作为第一相关状态参数Urelated
步骤D,构建用于同时编码上述第一相关坐标信息和上述第一相关状态参数的第一子量子线路Oracle Ob1,Oracle Ob1用于实现:
|i>|0>|0>→|i>|Urelated>|Xrelated>。
步骤E,根据Urelated和Xrelated获得bi(Xrelated,Urelated)。
步骤F,构建编码bi(Xrelated,Urelated)的第二子量子线路Oracle Ob2,Oracle Ob2用于实现|i>|Urelated>|Xrelated>|0>→|i>|Urelated>|Xrelated>|bi>。
步骤G,构建第三子量子线路Oracle Ob3,其中:上述第三子量子线路Oracle Ob3是上述第一子量子线路Oracle Ob1的转置共轭量子线路,Oracle Ob3用于实现|i>|Urelated>|Xrelated>|bi>→|i>|0>|0>|bi>。
步骤H,构建用于量子态转化的第四子量子线路Oracle Ob4,Oracle Ob4用于实现|i>|0>|0>|bi>→bi|i>。
步骤I,基于上述第一子量子线路Oracle Ob1、上述第二子量子线路Oracle Ob2、上述第三子量子线路Oracle Ob3和上述第四子量子线路Oracle Ob4构建Oracle Ob
需要说明的是,借助上述第一子量子线路Oracle Ob1、上述第二子量子线路OracleOb2、上述第三子量子线路Oracle Ob3和上述第四子量子线路Oracle Ob4构建Oracle Ob只是一种可能方式,一方面上述第一子量子线路Oracle Ob1、上述第二子量子线路Oracle Ob2、上述第三子量子线路Oracle Ob3和上述第四子量子线路Oracle Ob4并不遵守严格的顺序步骤限制,另一方面,每个子量子线路,其本质均是实现某种功能的一个量子态演化操作,即等效一个量子量子逻辑门。量子逻辑门在具体应用场景中,可以作用在用于描述具体应用场景的不同的量子态上和/或辅助量子比特位的预设量子态上,更直观的理解是,初始量子态的选择和设置不是唯一的,这是符合量子态叠加特性的,也是基于量子比特资源的占用和释放方便性,本领域技术人员可以按需进行设置的,示例性的,在|i>|Urelated>|Xrelated>|0>→|i>|Urelated>|Xrelated>|bi>的实现过程中,可以设置初始态是|i>|Urelated>|Xrelated>|0>之前的其它态,例如初始态是|i>|Urelated>|Xrelated>|1>。此时,为实现量子态演化操作的可逆操作,需要在通过额外量子线路进行量子比特资源的释放。
在本实施例中,优选设置初始态为|i>|Urelated>|Xrelated>|0>,一方面可以减少量子态初始态制备的复杂程度,另一方面,可以在实现该功能性的量子态演化操作的可逆操作的过程中,减少量子线路的在编码的复杂程度,整体可以起到简化量子线路的设置的效果。
同理,对上述其它功能Oracle的设置,也满足以上量子线路简单方便优化设置的需求。故可以理解的是,以上提供的步骤A到步骤I只是Oracle Ob的一种可能的构建方式。Oracle OA对网格单元的流体状态的作用可以为:
Figure BDA0002616833180000161
Ai′j′为系数矩阵的一个元素,
Figure BDA0002616833180000162
为异或符号,Z为自定义二进制数,可以为0,此时公式(11)可以表示为:OA|i>|j>=|i>|j>|Ai′j′>。
系数矩阵Ai′j′表示为:
Figure BDA0002616833180000163
i0、j0分别表示网格单元i、j的编号,即i、j的具体数值,每个网格单元有4个状态参数,i1、j1范围是0~3的整数(包括0及3),分别表示网格单元i0、j0的4种状态参数。
Figure BDA0002616833180000171
表示网格单元i0中的状态参数,,Uj0表示网格单元j0中的状态参数;
Figure BDA0002616833180000172
表示网格单元i0相关的位置坐标,Xj0表示网格单元j0相关的位置坐标。
这里函数
Figure BDA0002616833180000173
形式对于不同的i0,i1,j0,j1也是有规律的,它们可以分为2种类型:i0=j0和i0≠j0,对于不同的i1和j1,每种类型又包含16种情况,因此,可以有效地构建OA
为了计算
Figure BDA0002616833180000174
首先构造:
Figure BDA0002616833180000175
Figure BDA0002616833180000176
就是
Figure BDA0002616833180000177
的算术表达式,从而得到:
Figure BDA0002616833180000178
因为
Figure BDA0002616833180000179
Figure BDA00026168331800001710
为为了计算目标状态而构造出的辅助态,因此可以不必计算。
从而实现Oracle OA,其对应的量子线路如图8所示。所需的量子比特数是2m+CAqf,CAqf是用来计算Aij的辅助量子比特序列。量子逻辑门的复杂度为O(poly(m))。Aij的表达式比bi更复杂,所以CA≥Cb
下面对Oracle OA的一种可能的构建方式进行说明,在一种可能的实施方式中,上述基于上述第一量子线路、上述第二量子线路、上述量子随机存取存储器构建表示上述网格单元流体运动状态变化的线性系统方程的参数的第三量子线路,包括:
步骤1,针对An包括的系数矩阵Ai′j′,基于上述第一Oracle获取网格单元i的第j个相邻网格单元的网格数g(i,j),其中:g(i,j)为网格单元的状态参数和坐标信息的函数。
步骤2,基于网格单元i的第j个相邻网格单元的网格数g(i,j)构建Ai′j′的第i行第j个非零元素的列数f(i,j);其中,
Figure BDA0002616833180000181
其中:f(i,j)为网格单元的状态参数和坐标信息的函数,%表示取余运算,
Figure BDA00026168331800001811
表示向下取整符号。
此处量子线路的构造有一个技巧,例如图9所示。可以简单地对i和j的最高log(N)量子比特进行第一Oracle运算,然后,在|i>和|j>的整个空间中构造了Ol,其中Ol对网格单元状态的作用可以写为:Ol|i>|j>=|i>|f(i,j)>,在这个过程中所需的量子比特数是2m。
步骤3,根据f(i,j)获得Aij,其中:
Figure BDA0002616833180000182
i0、j0表示网格单元i、j的编号,每个网格单元有4个状态参数,i1、j1范围是0~3,分别表示网格单元i0、j0的4种状态参数;
Figure BDA0002616833180000183
表示网格单元i0中的状态参数,,Uj0表示网格单元j0中的状态参数;
Figure BDA0002616833180000184
表示网格单元i0相关的位置坐标,Xj0表示网格单元j0相关的位置坐标。
步骤4,针对i0和j0,利用上述量子随机存取存储器提取其对应的状态参数
Figure BDA0002616833180000185
和Uj0;并利用上述第二Oracle提取其对应的坐标信息
Figure BDA0002616833180000186
和Xj0
步骤5,构建用于同时编码上述坐标信息
Figure BDA0002616833180000187
和Xj0和上述状态参数
Figure BDA0002616833180000188
和Uj0的第一子量子线路Oracle OA1,Oracle OA1用于实现:
Figure BDA0002616833180000189
步骤6,构建编码
Figure BDA00026168331800001810
的第二子量子线路Oracle OA2,Oracle OA2用于实现:
Figure BDA0002616833180000191
步骤7,构建第三子量子线路Oracle OA3,其中:上述第三子量子线路Oracle OA3是上述第一子量子线路Oracle OA3的转置共轭量子线路,Oracle OA3用于实现
Figure BDA0002616833180000192
Figure BDA0002616833180000193
步骤8,基于上述第一子量子线路Oracle OA1、上述第二子量子线路Oracle OA2和上述第三子量子线路Oracle OA3构建Oracle OA
同上述借助第一子量子线路Oracle Ob1、上述第二子量子线路Oracle Ob2、上述第三子量子线路Oracle Ob3和上述第四子量子线路Oracle Ob4的构建Oracle Ob过程只是一种可能方式的论述,上述借助第一子量子线路Oracle OA1、上述第二子量子线路Oracle OA2和上述第三子量子线路Oracle OA3构建Oracle OA的过程也只是一种可能方式,在此不在赘述。
在一种可能的实施方式中,上述针对所有的网格单元,基于上述第三量子线路求解上述网格单元的线性系统方程,当上述网格单元的流体状态趋于稳定时,获得上述网格单元线性系统方程的状态参数表示的流体状态作为上述网格单元的目标状态,具体包括:
步骤一,针对所有的网格单元,基于上述第三量子线路表示的上述线性系统方程,从上述网格单元的状态参数的初始值开始,对上述线性系统方程进行迭代直至bn趋近于0时,获取此时的ΔUn
步骤二,根据获取的ΔUn对上述量子随机存取存储器中存储的Un进行更新,得到Un +1作为上述网格单元的目标状态。
将问题线性化后,转化为求解线性方程,该过程的目标是输出一个目标状态:
Figure BDA0002616833180000201
这样
Figure BDA0002616833180000202
Figure BDA0002616833180000203
为计算得到的目标状态,|x>为目标状态的真值,|x>的振幅满足
Figure BDA0002616833180000204
其中,∈为预设的精度。
此处可以使用安德鲁·柴尔德(Andrew Childs)的Chebyshev方法的量子线性求解器(quantum linear solver)作为计算的加速器,该量子线性求解器的整体量子线路可以如图10所示。其中,t=O(log(j0))是控制比特
Figure BDA0002616833180000205
的量子比特数,α为
Figure BDA0002616833180000206
的归一化常数。总共的量子比特数是2(m+1)+t+CAqf。通过量子线性求解器得到
Figure BDA0002616833180000207
Figure BDA0002616833180000208
转化为数字形式:
Figure BDA0002616833180000209
根据MU生成以下状态:
Figure BDA00026168331800002010
MU为量子随机存储器中存储的所有网格单元的流体状态。
接着在
Figure BDA00026168331800002011
Figure BDA00026168331800002012
上执行量子加法器,得到:
Figure BDA00026168331800002013
Figure BDA00026168331800002014
解纠缠
Figure BDA00026168331800002015
得到
Figure BDA00026168331800002016
其量子线路如图11所示。
从而得到
Figure BDA00026168331800002017
并用
Figure BDA00026168331800002018
更新MU。在此过程中,需要m+2qf个量子比特。用到两次MU,一次QADC和一次量子加法器。
从量子随机存取存储器中获取演化得到的态为数字编码态
Figure BDA00026168331800002019
即网格单元的目标状态。
在一种可能的实施方式中,上述ΔUn为信息编码在量子态振幅上的模拟编码态数据;上述根据ΔUn对上述量子随机存取存储器中存储的Un进行更新,得到Un+1作为上述网格单元的目标状态,具体包括:
步骤1,构建第一数据转化量子线路Oracle Oconvert1,Oracle Oconvert1用于实现:
ΔUn|i>|0>→|i>|ΔUn>。
步骤2,从上述量子随机存取存储器中提取Un
步骤3,构建量子加法器线路Oracle Oadd,Oracle Oadd用于实现:
|i>|Un>|ΔU>→|i>|Un>|Un+1>。
步骤4,构建第二数据转化量子线路Oracle Oconvert2,Oracle Oconvert2用于实现:|i>|Un>|Un+1>→|i>|0>|Un+1>。
步骤5,获得由依次设置的上述第一数据转化量子线路Oracle Oconvert1、上述量子随机存取存储器、上述量子加法器线路Oracle Oadd和上述第二数据转化量子线路OracleOconvert2的量子线路编码的由量子态本征态编码的|i>|0>|Un+1>作为上述线性系统方程所包含的状态参数目标值。
步骤6,根据上述状态参数目标值得到上述网格单元的目标状态。
需要说明的是,同上述借助第一子量子线路Oracle Ob1、上述第二子量子线路Oracle Ob2、上述第三子量子线路Oracle Ob3和上述第四子量子线路Oracle Ob4的构建Oracle Ob过程只是一种可能方式的论述,上述借助上述第一数据转化量子线路OracleOconvert1、上述量子随机存取存储器、上述量子加法器线路Oracle Oadd和上述第二数据转化量子线路Oracle Oconvert2得到量子线路的过程也只是一种可能方式,在此不在赘述。在一种可能的实施方式中,上述根据上述状态参数目标值得到上述网格单元的目标状态,具体包括:
步骤a,构建第三数据转化量子线路Oracle Oconvert3,Oracle Oconvert3用于实现:|i>|0>|Un+1>→Un+1|i>。
需要说明的是,Oracle Oconvert3用于实现数字编码态到模拟编码态的转变,但是所转变的量子态并不是唯一的,示例性的:Oracle Oconvert3用于实现:|i>|1>|Un+1>→Un+1|i>|1>,此时对应实现Oracle Oconvert3的量子线路将需要多一位的辅助量子比特编码|1>态,不满足量子线路设置时需要尽量简化的需求。
步骤b,利用上述第三数据转化量子线路Orracle Oconvert3将上述状态参数目标值转换为由量子态振幅编码的Un+1,获得各上述网格单元的目标状态。
线性系统AnΔUn=bn能够通过量子算法来解。然而,输出结果是相关的态,而不是真实解的向量形式。不同于向量
Figure BDA0002616833180000221
量子态|ΔU>已经被归一化:
Figure BDA0002616833180000222
对于不同的A和b,有不同的归一化常数;通常第n次迭代的常数cn不同于第n+1次迭代的常数cn+1
Figure BDA0002616833180000223
等于
Figure BDA0002616833180000224
也就是说执行像
Figure BDA0002616833180000225
的操作是必须的。还需要实现|ΔUi>加上|ΔUj>。假设存在充分大的常数L,构建一个新的态|ΔU′>:
Figure BDA0002616833180000226
其中,
Figure BDA0002616833180000227
进而得到方程的一个新集合:
Figure BDA0002616833180000228
当Δt充分小时,线性化方程(19)中的第二个方程,用Δu′i n-1Δu′i n提换Δu′i 2。从而得到一个新的线性系统
Figure BDA0002616833180000229
其中,
Figure BDA0002616833180000231
对于态|ΔU′>,归一化常数总是1。然后,就能够自由地完成数学操作(用|ΔU′>的前N部分)。当N充分大时,A′也是一个稀疏矩阵,因此能成功构建它。
可选的,本申请实施例基于量子算法的计算流体动力学模拟方法可以应用到SU2中,通过量子子程序替换某些模块,在执行过程中寻求量子加速。在SU2中替换的模块为CfluidIteration(计算流体迭代)类里面使用到的数据结构和方法。在数据结构中,它们被量子存储器或量子状态的集合所取代,这些量子状态包括关于网格设置和线性方程中涉及的所有向量的信息。对于时间积分,线性解算器被量子解算器所取代。在预处理和后处理阶段,它们完全被量子算法模块所取代,形成量子线性求解器所必需的Oracle。此外,并没有改变问题模型、网格设置和有限体积法,包括时间积分的迭代过程。
量子算法的时间复杂度与logN成正比,其中N表示网格中的单元数或网格数。量子比特的数量也与logN成比例。此外,网格数量、复杂性和量子比特的数量都受到精度∈的影响;条件数κ和稀疏性s也会影响时间复杂性,可选的,因为精度∈、条件数κ及稀疏性s的变化不大,所以在计算中可以认为他们是不变的。
SU2是一个开源的流体动力学解决方案,本申请实施例中利用量子子过程代替经典求解器中的一些模块,形成一个量子-经典混合模型,最后演示量子加速,其中,经典模型例如图12a所示,量子-经典混合模型例如图12b所示。
针对SU2_CFD模块,程序从读取配置文件的主函数(main)开始。然后通过计算流体驱动器(CfluidDriver)实例来控制整个计算过程,即执行计算流体迭代。时间积分在一个由CfluidIteration控制的for循环中执行。在每个迭代步骤中,调用MultiGrid_Cycle(多重网格循环)方法和预处理方法。MultiGrid_Cycle包含计算FVM问题的MultiGrid(多重网格)算法。预处理连接上一个迭代的结果和下一个迭代的输入。在MultiGrid_Cycle中,预处理子进程读取网格的数据并生成Jacobi(雅克比)矩阵,然后将其发送到CeulerSolver(流体控制方程处理器)并通过经典线性求解器求解线性方程。
CmultiGridIntegration(多重网格积分)的预处理模块包括预处理生成的线性问题。这部分被一个对应的量子算法所替代,并调用包含网格信息数据的Oracle。与经典过程相似,最终目标是生成下一步所需的量子数据:量子线性求解器。
模块MultiGrid_Cycle调用CMultiGridIntegration中的TimeIntegration(时间积分)模块,取而代之的是量子线性求解器,它执行任务的速度要比经典解算器快得多。
CFluidIteration中的预处理处理MultiGrid_Cycle之后的结果,然后为下一次迭代生成数据,并显示结果。这个模块被一个量子后处理所取代,其目的是向向量中添加变化并准备下一次迭代。
量子计算过程可以如图13所示,其中,图12a中标识101的部分被替换为QuantumPreprocessing(量子预处理),标识102的部分被替换为Quantum Linear Solver(量子线性求解器),标识103的部分被替换为Quantum PostProcessinig(量子后处理)。每个迭代步骤中的量子过程都会执行一步时间积分Δt,输入
Figure BDA0002616833180000241
输出
Figure BDA0002616833180000242
有如下关系:
Figure BDA0002616833180000243
Figure BDA0002616833180000244
由此可以得到量子线性解算器所需的量子Oracles。这里在第i+1次迭代中,Ai+1表示系数矩阵,bi+1表示残差量,
Figure BDA0002616833180000245
表示网格单元的坐标信息的量子Oracle。然后用量子计算机求解如下方程:
Figure BDA0002616833180000246
然后,通过在最后一步
Figure BDA0002616833180000247
中增加变化量,这个迭代步骤就结束了:
Figure BDA0002616833180000248
其具体实现过程可以参见上述实施例中的相关内容。在模型中QRAM是演示量子加速的必要条件,在三个过程中使用QRAM:一是存储关于网格单元的信息,二是存储求解器的区间结果,三是帮助从执行非酉转换的概率错误中恢复。该算法整个过程中量子比特、查询和量子逻辑门的复杂度。所需的量子比特数是:
Figure BDA0002616833180000251
网格单元数量与量子比特数的关系如图14所示,QRAM和Oracle查询的复杂性为:
Figure BDA0002616833180000252
量子逻辑门的复杂度为:
Figure BDA0002616833180000253
将SU2和量子编程框架相结合,编写了求解CFD问题的量子-经典混合求解程序,可以使用QPanda量子编程框架,QPanda及其OriginIR(量子中间表示)支持任何量子Qracle的定义和模拟。
通过编程实现本申请的量子线性求解器部分,该部分采用生成的线性方程,以SU2为输入,并输出量子线路编码方程式。该程序结合了SU2模块和量子编程框架QPanda。使用量子Oracles进行处理,QPanda支持定义和模拟量子Oracles。量子Oracle是一个黑盒表示某种量子态的转变。量子Oracles的一个典型例子是量子函数:O|x>|0>=|x>|f(x)>,此处计算f(x)用第一个量子寄存器作为输入,第二个量子寄存器作为输出。另外一个例子就是QRAM可以被视为一种Oracle。很多量子算法用到Oracle,许多量子算法都使用Oracle,但他们不关心Oracle的实现-它可以分解成量子门,也可以实现QRAM。在QPanda中,可以使用“Oracle”函数来定义。Oracle被认为是具有用户提供的名称。
网格生成模块除了可以加速时间积分过程外,还可以集成在量子计算机中。在上述算法中,将网格信息视为Oracle。一种直接的方法是预先生成网格文件,将经典数据编码到量子计算机中,也可以用量子算法代替网格生成程序,而不是读取预先生成的网格文件。求解器需要一个有效的矩阵A和常数
Figure BDA0002616833180000254
的量子表示。在这个问题中,偏微分方程组需要被线性化,以得到A和
Figure BDA0002616833180000255
的表达式。
为了提升本申请实施例中的加速效果,量子计算机需要满足以下要求:容错,本申请实施例的方法加快了每个迭代步骤的计算速度,然而,随着网格数量的增加,迭代的次数会大大增加。过程中不执行量子测量来提取每一步的结果。而在整个计算过程中保留量子特征(纠缠、叠加)。该方法的这一特点意味着量子态的存储不受退相干和噪声的影响,具有容错能力的量子计算机能够增加计算结果的准确性。
量子-经典混合算法面临着经典与量子之间的数据转换问题。问题是在∑|addr>|data>和经典内存之间执行这种转换的时间复杂度是多少。从量子到经典的一个简单方法是重复测量,或者是经常提到的“quantum tomography”。一般情况下从经典到量子的方式是通过一组受控的旋转操作在量子线路中对状态准备阶段进行编码。这些方法很简单,但是会消耗O(2n)时间。如果不能有效地实现量子经典转换,许多量子算法将永远无法实现次线性时间复杂度,因为输入和输出至少会消耗线性时间。为了解决这个问题,需要使用量子随机存取存储器。QRAM是一种可以有效地存储和加载量子信息的结构,QRAM在叠加中操作地址。它使数字形式的量子数据能够与时间复杂度O(logN)并行加载和存储。QRAM是经典随机存取存储器的量子模拟。当地址处于叠加状态时,QRAM使用户能够有效地获取QRAM中的数据。
QRAM用于实现:∑i|i>|0>=∑i|i>|Ui>。
其中,Ui表示第i个网格单元的量子态编码的经典数据。
QRAM模型可以使用一组光学设备实现,需要O(logN)操作(N表示内存的长度)。QRAM模型也通过分解量子线路来实现。
在申请实施例中,构建了混合量子-经典方法的量子模块替换SU2 CFD解算器中的某些模块,能够是实现通过计算流体动力学模拟流体运动。原始的经典求解器使用FVM以及迭代时间积分来获得流体方程的解。通过多重网格方法,经典求解器在每个迭代步骤之间具有O(n)复杂度;这里的n表示网格数。本申请实施例中使用量子线性求解器来加速线性方程组的求解,替换CFluidIteration过程中的核心过程。一个是MultiGrid_Cycle,包括预处理和方程求解,已被量子过程完全取代。另一个是预处理,它在两个迭代步骤之间链接数据,并且也被量子过程所取代。经典数据自然地被量子数据取代。输入的网格信息应放入预处理器中经常使用的QRAM中。不再需要计算过程中的大多数内部数据,因为它们被保存在量子寄存器中。整个问题的配置数据被视为常数,用于直接生成量子程序。本申请实施例的方案中时间复杂度O(log(N)),N为网格数,计算速度比经典对数快得多。
在申请实施例中,在不更改问题模型的情况下,实现了该问题的指数量子加速。其次,充分展示了SU2程序,找到必须做的适当替换,然后借助QPanda量子编程框架实现我们的算法,最终获得了数值结果。并且不忽略用于数据转换的量子-经典接口问题。
以下对本申请实施例中的量子线路构造中使用的量子线路模块进行介绍。
1.首先对QDAC(Quantum Digital-to-Analog Converter,量子数模转换器)及QADC(Quantum Analog-to-Digital Converter,量子模数转换器)这两个量子线路模块进行介绍。
QDAC的目的是改变由0到1之间的数字编码的量子态到量子态对应的振幅上。因此,QDAC操作由定义:
Figure BDA0002616833180000271
QDAC算法的主要思想是通过旋转辅助位来在每个基态前面转换系数,以下从三个方面进行解释。
首先,解释成功的可能性。QDAC算法只有一定的成功概率和一定的操作失败概率。如果失败,请从头开始运行算法。经过计算,一次成功的概率为
Figure BDA0002616833180000272
其次,QDAC算法的量子线路如图15所示。在测量时如果是|0>态,则继续;如果是|1>态,则从头开始做一遍。n表示编码数据所需的量子比特数,m表示所需精度的位数,其中只需要计算高位。UD门表示为:UD({dj})|j>|0>=|j>|dj>。UAC门代表arccos函数,受控Ry门表示一个比特的受控旋转。如果不考虑n个H门,则一共需要O(m)个门。
然后,解释归一化问题。最终系数需要归一化。经过计算,归一化系数C是
Figure BDA0002616833180000273
QADC算法包括绝对值QADC,实部QADC和虚部QADC,这些算法的主要思想涉及Hadamard测试和相位估计。本申请中主要运用实部QADC。一个m比特实部QADC操作被定义为:
Figure BDA0002616833180000274
其中,
Figure BDA0002616833180000275
表示m比特串
Figure BDA0002616833180000276
使得
Figure BDA0002616833180000277
最接近于ck的实部。实部QADC算法的量子线路如图16a所示。相位估计的量子线路如图16b虚线框中所示。实部QADC的量子线路Q部分如图16c虚线框中所示。m和n分别代表所需精度的位数和编码数据所需的位数。一个模拟编码酉变换UA({cj})被定义为:
Figure BDA0002616833180000281
这里
Figure BDA0002616833180000282
是经典数据,并且满足归一化条件
Figure BDA0002616833180000283
门G被下式给出:
Figure BDA00026168331800002812
其中,Q是图16c中显示的门,S0=I-2(|0><0|)data,B是一个有条件的相位翻转门,ZB是仅作用在量子比特B上的泡利Z门。PE表示相位估计,IQFT是指逆量子傅里叶变换。
在相位估计后,Q的相位θk可以被编码到register0量子比特上,然后使用量子算术运算Uf(θ)去计算xk=f(θk)=2sin2πθk-1。最后,未计算的数据,B,register0量子比特,数字编码状态
Figure BDA0002616833180000284
被得到。在实部QADC量子线路中,需要2(m+n)+1个量子比特。由于门S0可被分解为
Figure BDA0002616833180000285
个单门和双门。因此,单和双量子位门的数量是
Figure BDA0002616833180000286
O(1/∈)个受控UA门被使用。
2.通过操作线性组合实现矩阵逆。
任何简单的函数都可以线性近似为其他函数的线性组合,可以通过Chebyshev多项式近似矩阵的逆函数。具体如下:
Figure BDA0002616833180000287
Figure BDA0002616833180000288
Figure BDA0002616833180000289
还有
Figure BDA00026168331800002810
本申请中需要求逆矩阵函数f(x)=1/x,它满足O(‖A-1-F‖)=∈。线性组合为:
Figure BDA00026168331800002811
此处
Figure BDA0002616833180000291
g(x)是2∈,
在Dκ:=(-1,-1/κ)∪(1/κ,1)。
Figure BDA0002616833180000292
是第一类Chebyshev多项式。
量子游走:为了实现Chebyshev多项式,需要在量子行走框架中进行。
因为量子游走被执行在空间
Figure BDA0002616833180000293
中的一些态
Figure BDA0002616833180000294
上,定义一个映射
Figure BDA0002616833180000295
从从
Figure BDA0002616833180000296
Figure BDA0002616833180000297
Figure BDA0002616833180000298
和行走算子:
Figure BDA00026168331800002919
算子S执行
Figure BDA0002616833180000299
中的乘积态的翻转操作。于是有:
Figure BDA00026168331800002910
Figure BDA00026168331800002911
是第一类Chebyshev多项式。
量子线路实施:用两个Oracles来提取矩阵A的信息:Ol和OA,写作:
Figure BDA00026168331800002912
第一个目标是去编码态
Figure BDA00026168331800002913
可以准备O(j0)长度的量子线路。
接下来,构建T算子,算子T可以如下合成。对于向量|j>,整个线路需要
Figure BDA00026168331800002914
个量子比特,同样,
Figure BDA00026168331800002915
表示其乘积状态和一堆辅助量子位用作控制旋转。
Figure BDA00026168331800002916
Figure BDA00026168331800002917
Figure BDA00026168331800002918
把Ajl看作一系列控制比特,定义一个控制旋转算子:
Figure BDA0002616833180000301
Figure BDA0002616833180000302
算子M被一个Oracle Oamp合成为:
Figure BDA0002616833180000303
目的是准备一个对目标态的逼近。
Figure BDA0002616833180000304
关键的贡献是作为控制旋转算子的M。
|α>|0〉:=|α>(sin(θ)|0>+cos(θ)|1>)
再次运用OH把|Hjl>变回|0>。
下面,例如图17a-17c所示,说明W的量子线路。因为
Figure BDA0002616833180000305
S可以由一群交换操作构造出来,剩下来的就是
Figure BDA0002616833180000306
Figure BDA0002616833180000307
所以
Figure BDA0002616833180000308
算子T在量子线路中是酉算子,它把|j>|0>变为|ψj>。它不同于维数为4N2×2N的T:=∑j∈Nj><j|。为了区别,用Tqc来定义该量子线路:
Figure BDA0002616833180000309
在量子线路中,H表示H门,OH、M表示不同功能的Oracle,
Figure BDA00026168331800003011
表示转置共轭,T表示H门及Oracle组合的整体功能模块T,T模块的功能即为将|j>变换为|ψj>。并且,获得的输入该T模块的矩阵为N阶矩阵,上方|j>处
Figure BDA00026168331800003010
中的N表示矩阵行数,下方
Figure BDA0002616833180000311
中的N表示矩阵列数,其余表示同上。构建的T模块在量子线路中可等效于量子逻辑门,其矩阵形式为:
Figure BDA0002616833180000312
其中,
Figure BDA0002616833180000313
为量子态左矢。
具体的,利用H门,构建叠加态:
Figure BDA0002616833180000314
其中,j为目标取值,表示矩阵的第j行;
Figure BDA0002616833180000315
表示张量积或张乘;d为第j行的非0元素总数;l为非0元素在第j行所有非0中元素的序号,表示第l个非0元素,|l>对应的量子比特位即为第一比特位。
OH实现变换:
Figure BDA0002616833180000316
其中,Hjl′为矩阵第j行第l′列的非0元素数值,|Hjl′>对应的量子比特位即为第二比特位。需要说明的是,如果Hjl′为复数,可以将实部和虚部编码到第二比特位上,即|Hjl′>=|real>|imag>,real表示实部,image表示虚部;如果Hjl′写成欧拉形式re,则可将r和θ的信息编码到第二比特位,即|Hjl′>=|r>|θ>。
M实现变换:
Figure BDA0002616833180000317
最后,再调用一次OH进行转置共轭操作,将编码Hjl′的量子比特回收,进而输出|ψj>。其中,Hjl′ *为矩阵H第j行第l′列的非0元素数值的共轭,Hmax为矩阵H中绝对值最大的元素数值。
本申请实施例还提供了一种基于量子算法的计算流体动力学模拟装置,参见图18,该装置包括:
网格量子线路构建模块201,用于在利用有限体积法解析计算流体动力学的过程中,针对流体运动的离散化数值网格中各网格单元,构建表示上述网格单元坐标信息的第一量子线路、表示上述网格单元的状态参数的第二量子线路,其中,上述网格单元的状态参数存储在量子随机存取存储器中,上述量子随机存取存储器可操作处于量子叠加态的地址和数据;
参数量子线路构建模块202,用于基于上述第一量子线路、上述第二量子线路、上述量子随机存取存储器构建表示上述网格单元流体状态变化的线性系统方程的参数的第三量子线路;
目标状态获取模块203,用于针对所有的网格单元,基于上述第三量子线路求解上述网格单元的线性系统方程,当上述网格单元的流体状态趋于稳定时,获得上述网格单元线性系统方程的状态参数表示的流体状态作为上述网格单元的目标状态。
在一种可能的实施方式中,所述第一量子线路包括第一Oracle及第二Oracle;所述第一Oracle用于提取指定网格单元的相邻网格数,所述第二Oracle用于提取网格单元的坐标信息。
在一种可能的实施方式中,所述第二量子线路,包括:通过量子态的振幅提取所述网格单元的状态参数的第二模拟态量子线路,和通过量子态的本征态提取所述网格单元的状态参数的第二数字态量子线路。
在一种可能的实施方式中,所述量子随机存取存储器内存储有的地址和数据的预设映射关系,所述地址用第一叠加态编码,所述第一叠加态中的每个本征态对应一个地址信息,所述数据用第二叠加态编码,所述第二叠加态中的每个本征态对应一个数据信息。
在一种可能的实施方式中,所述线性系统方程为:AnΔUn=bn;其中:ΔUn=Un+1-Un;其中:n表示记录流体状态的第n个时刻;Un表示第n个时刻网格单元的流体状态,Un={Ui},i=0,1,2,…,N-1,Ui表示网格单元i的流体状态;An表示第n个时刻的系数矩阵,且
Figure BDA0002616833180000321
表示An为与
Figure BDA0002616833180000322
相关的函数,bn表示第n个时刻的残差量,且
Figure BDA0002616833180000323
表示bn为与
Figure BDA0002616833180000324
相关的函数,N为网格单元的总数量,
Figure BDA0002616833180000325
为所述网格单元的坐标信息;
所述线性系统方程的参数包括当前时刻的系数矩阵An和当前时刻的残差量bn;所述第三量子线路包括对应所述系数矩阵An的Oracle OA,和对应所述残差量bn的Oracle Ob;所述Oracle OA用于提取当前n时刻系数矩阵An的元素Ai′j′,An={Ai′j′},i′=4i+i1,j′=4j+j1,i,j=0,1,2,…,N-1;i1,j1是0到3的整数,对网格单元i的第j个相邻网格单元,所述Oracle OA的作用为OA|i>|j>=|i>|j>|Ai′j′>;所述Oracle Ob用于提取当前n时刻的残差量bn的元素bi,其中,bi表示网格单元i的残差量,bn={bi},i=0,1,2,…,N-1;对网格单元i,所述Oracle Ob的作用为Ob|0>=∑ibi|i>。
在一种可能的实施方式中,所述参数量子线路构建模块202,具体用于:
针对bn包括的对应于网格单元i的残差量bi,获取网格单元i的相关网格单元;
基于所述第二Oracle获取所述相关网格单元的坐标信息作为第一相关坐标信息Xrelated
基于所述第二量子线路和所述量子随机存取存储器,获取相关网格单元对应的状态参数信息作为第一相关状态参数Urelated
构建用于同时编码所述第一相关坐标信息和所述第一相关状态参数的第一子量子线路Oracle Ob1,Oracle Ob1用于实现:
|i>|0>|0>→|i>|Urelated>|Xrelated>;
根据Urelated和Xrelated获得bi(Xrelated,Urelated);
构建编码bi(Xrelated,Urelated)的第二子量子线路Oracle Ob2,Oracle Ob2用于实现|i>|Urelated>|Xrelated>|0>→|i>|Urelated>|Xrelated>|bi>;
构建第三子量子线路Oracle Ob3,其中:所述第三子量子线路Oracle Ob3是所述第一子量子线路Oracle Ob1的转置共轭量子线路,Oracle Ob3用于实现|i>|Urelated>|Xrelated>|bi>→|i>|0>|0>|bi>;
构建用于量子态转化的第四子量子线路Oracle Ob4,Oracle Ob4用于实现|i>|0>|0>|bi>→bi|i>;
基于所述第一子量子线路Oracle Ob1、所述第二子量子线路Oracle Ob2、所述第三子量子线路Oracle Ob3和所述第四子量子线路Oracle Ob4构建
Oracle Ob
在一种可能的实施方式中,所述参数量子线路构建模块202,具体用于:
针对An包括的系数矩阵Ai′j′,基于所述第一Oracle获取网格单元i的第j个相邻网格单元的网格数g(i,j),其中:g(i,j)为网格单元的状态参数和坐标信息的函数;
基于网格单元i的第j个相邻网格单元的网格数g(i,j)构建Ai′j′的第i行第j个非零元素的列数f(i,j);其中,
Figure BDA0002616833180000341
其中:f(i,j)为网格单元的状态参数和坐标信息的函数,%表示取余运算;
Figure BDA0002616833180000342
表示向下取整符号;
根据f(i,j)获得Ai′j′,其中:
Figure BDA0002616833180000343
i0、j0分别表示网格单元i、j的编号,每个网格单元有4个状态参数,i1、j1范围是0~3,分别表示网格单元i0、j0的4种状态参数;
Figure BDA0002616833180000344
表示网格单元i0中的状态参数,,Uj0表示网格单元j0中的状态参数;
Figure BDA0002616833180000345
表示网格单元i0相关的位置坐标,Xj0表示网格单元j0相关的位置坐标;
针对i0和j0,利用所述量子随机存取存储器提取其对应的状态参数
Figure BDA0002616833180000351
和Uj0;并利用上述第二Oracle提取其对应的坐标信息
Figure BDA0002616833180000352
和Xj0
构建用于同时编码所述坐标信息
Figure BDA0002616833180000353
和Xj0和所述状态参数
Figure BDA0002616833180000354
和Uj0的第一子量子线路Oracle OA1,Oracle OA1用于实现:
Figure BDA0002616833180000355
构建编码
Figure BDA0002616833180000356
的第二子量子线路Oracle OA2,Oracle OA2用于实现:
Figure BDA0002616833180000357
构建第三子量子线路Oracle OA3,其中:所述第三子量子线路Oracle OA3是所述第一子量子线路Oracle OA3的转置共轭量子线路,Oracle OA3用于实现:
Figure BDA0002616833180000358
基于所述第一子量子线路Oracle OA1、所述第二子量子线路Oracle OA2和所述第三子量子线路Oracle OA3构建Oracle OA
在一种可能的实施方式中,所述目标状态获取模块203,包括:
迭代计算子模块,用于针对所有的网格单元,基于所述第三量子线路表示的所述线性系统方程,从所述网格单元的状态参数的初始值开始,对所述线性系统方程进行迭代直至bn趋近于0时,获取此时的ΔUn
目标状态计算子模块,用于根据获取的ΔUn对所述量子随机存取存储器中存储的Un进行更新,得到Un+1作为所述网格单元的目标状态。
在一种可能的实施方式中,所述ΔUn为信息编码在量子态振幅上的模拟编码态数据;所述目标状态计算子模块,包括:
第一构建单元,用于构建第一数据转化量子线路Oracle Oconvert1,Oracle Oconvert1用于实现:
ΔUn|i>|0>→|i>|ΔUn>;
状态提取单元,用于从所述量子随机存取存储器中提取Un
第二构建单元,用于构建量子加法器线路Oracle Oadd,Oracle Oadd用于实现:
|i>|Un>|ΔU>→|i>|Un>|Un+1>;
第三构建单元,用于构建第二数据转化量子线路Oracle Oconvert2,Oracle Oconvert2用于实现:|i>|Un>|Un+1>→|i>|0>|Un+1>;
目标值确定单元,用于获得由依次设置的所述第一数据转化量子线路OracleOconvert1、所述量子随机存取存储器、所述量子加法器线路Oracle Oadd和所述第二数据转化量子线路Oracle Oconvert2的量子线路编码的由量子态本征态编码的|i>|0>|Un+1>作为所述线性系统方程所包含的状态参数目标值;
目标状态确定单元,用于根据所述状态参数目标值得到所述网格单元的目标状态。
在一种可能的实施方式中,所述目标状态确定单元,具体用于:构建第三数据转化量子线路Oracle Oconvert3,Oracle Oconvert3用于实现:|i>|0>|Un+1>→Un+1|i>;利用所述第三数据转化量子线路Oracle Oconvert3将所述状态参数目标值转换为由量子态振幅编码的Un+1,获得各所述网格单元的目标状态。
本申请实施例还提供了一种量子计算机设备,包括:量子线路、量子效应器及量子随机存取存储器,该量子计算机设备在运行时实现上述任一所述的基于量子算法的计算流体动力学模拟方法。
本申请实施例还提供了一种计算机可读存储介质,该计算机可读存储介质内存储有计算机程序,上述计算机程序被处理器执行时实现上述任一所述的基于量子算法的计算流体动力学模拟方法。
需要说明的是,在本文中,各个可选方案中的技术特征只要不矛盾均可组合来形成方案,这些方案均在本申请公开的范围内。诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来,而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且,术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含,从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素,而且还包括没有明确列出的其他要素,或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下,由语句“包括一个……”限定的要素,并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。
本说明书中的各个实施例均采用相关的方式描述,各个实施例之间相同相似的部分互相参见即可,每个实施例重点说明的都是与其他实施例的不同之处。尤其,对于装置、量子计算机设备及存储介质的实施例而言,由于其基本相似于方法实施例,所以描述的比较简单,相关之处参见方法实施例的部分说明即可。
以上所述仅为本申请的较佳实施例而已,并非用于限定本申请的保护范围。凡在本申请的精神和原则之内所作的任何修改、等同替换、改进等,均包含在本申请的保护范围内。

Claims (13)

1.一种基于量子算法的计算流体动力学模拟方法,其特征在于,所述方法包括:
在利用有限体积法解析计算流体动力学的过程中,针对流体运动的离散化数值网格中各网格单元,构建表示所述网格单元坐标信息的第一量子线路、表示所述网格单元的状态参数的第二量子线路,其中,所述网格单元的状态参数存储在量子随机存取存储器中,所述量子随机存取存储器可操作处于量子叠加态的地址和数据;
基于所述第一量子线路、所述第二量子线路、所述量子随机存取存储器构建表示所述网格单元流体状态变化的线性系统方程的参数的第三量子线路;
针对所有的网格单元,基于所述第三量子线路求解所述网格单元的线性系统方程,当所述网格单元的流体状态趋于稳定时,获得所述网格单元线性系统方程的状态参数表示的流体状态作为所述网格单元的目标状态。
2.根据权利要求1所述的方法,其特征在于,所述第一量子线路包括第一Oracle及第二Oracle;
所述第一Oracle用于提取指定网格单元的相邻网格数,所述第二Oracle用于提取网格单元的坐标信息。
3.根据权利要求1所述的方法,其特征在于,所述第二量子线路,包括:
通过量子态的振幅提取所述网格单元的状态参数的第二模拟态量子线路,和通过量子态的本征态提取所述网格单元的状态参数的第二数字态量子线路。
4.根据权利要求1所述的方法,其特征在于,所述量子随机存取存储器内存储有的地址和数据的预设映射关系,所述地址用第一叠加态编码,所述第一叠加态中的每个本征态对应一个地址信息,所述数据用第二叠加态编码,所述第二叠加态中的每个本征态对应一个数据信息。
5.根据权利要求2所述的方法,其特征在于,所述线性系统方程为:AnΔUn=bn;其中:ΔUn=Un+1-Un
其中:n表示记录流体状态的第n个时刻;Un表示第n个时刻网格单元的流体状态,Un={Ui},i=0,1,2,...,N-1,Ui表示网格单元i的流体状态;An表示第n个时刻的系数矩阵,且
Figure FDA0002616833170000021
表示An为与Un
Figure FDA0002616833170000022
相关的函数,bn表示第n个时刻的残差量,且
Figure FDA0002616833170000023
表示bn为与Un
Figure FDA0002616833170000024
相关的函数,N为网格单元的总数量,
Figure FDA0002616833170000025
为所述网格单元的坐标信息;
所述线性系统方程的参数包括当前时刻的系数矩阵An和当前时刻的残差量bn
所述第三量子线路包括对应所述系数矩阵An的Oracle OA,和对应所述残差量bn的Oracle Ob
所述Oracle OA用于提取当前n时刻系数矩阵An的元素Ai′j′,An={Ai′j′},i′=4i+i1,j′=4j+j1,i,j=0,1,2,...,N-1;i1,j1是0到3的整数,对网格单元i的第j个相邻网格单元,所述Oracle OA的作用为OA|i〉|j〉=|i〉|j>|Ai′j′〉;
所述Oracle Ob用于提取当前n时刻的残差量bn的元素bi,其中,bi表示网格单元i的残差量,bn={bi},i=0,1,2,...,N-1;对网格单元i,所述Oracle Ob的作用为Ob|0>=∑ibi|i>。
6.根据权利要求5所述的方法,其特征在于,所述基于所述第一量子线路、所述第二量子线路、所述量子随机存取存储器构建表示所述网格单元流体状态变化的线性系统方程的参数的第三量子线路,包括:
针对bn包括的对应于网格单元i的残差量bi,获取网格单元i的相关网格单元;
基于所述第二Oracle获取所述相关网格单元的坐标信息作为第一相关坐标信息Xrelated
基于所述第二量子线路和所述量子随机存取存储器,获取相关网格单元对应的状态参数信息作为第一相关状态参数Urelated
构建用于同时编码所述第一相关坐标信息和所述第一相关状态参数的第一子量子线路Oracle Ob1,Oracle Ob1用于实现:
|i>|0>|0>→|i>|Urelated>|Xrelated>;
根据Urelated和Xrelated获得bi(Xrelated,Urelated);
构建编码bi(Xrelated,Urelated)的第二子量子线路Oracle Ob2,Oracle Ob2用于实现|i>|Urelated>|Xrelated>|0>→|i>|Urelated>|Xrelated>|bi>;
构建第三子量子线路Oracle Ob3,其中:所述第三子量子线路Oracle Ob3是所述第一子量子线路Oracle Ob1的转置共轭量子线路,Oracle Ob3用于实现|i>|Urelated〉|Xrelated>|bi>→|i>|0>|0>|bi>;
构建用于量子态转化的第四子量子线路Oracle Ob4,Oracle Ob4用于实现|i〉|0〉|0〉|bi〉→bi|i>;
基于所述第一子量子线路Oracle Ob1、所述第二子量子线路Oracle Ob2、所述第三子量子线路Oracle Ob3和所述第四子量子线路Oracle Ob4构建Oracle Ob
7.根据权利要求5或6所述的方法,其特征在于,所述基于所述第一量子线路、所述第二量子线路、所述量子随机存取存储器构建表示所述网格单元流体运动状态变化的线性系统方程的参数的第三量子线路,包括:
针对An包括的系数矩阵Ai′j′,基于所述第一Oracle获取网格单元i的第j个相邻网格单元的网格数g(i,j),其中:g(i,j)为网格单元的状态参数和坐标信息的函数;
基于网格单元i的第j个相邻网格单元的网格数g(i,j)构建Ai′j′的第i行第j个非零元素的列数f(i,j);其中,
Figure FDA0002616833170000041
其中:f(i,j)为网格单元的状态参数和坐标信息的函数,%表示取余运算;
Figure FDA0002616833170000048
表示向下取整符号;
根据f(i,j)获得Ai′j′,其中:
Figure FDA0002616833170000042
i0、j0分别表示网格单元i、j的编号,每个网格单元有4个状态参数,i1、j1范围是0~3,分别表示网格单元i0、j0的4种状态参数;
Figure FDA0002616833170000049
表示网格单元i0中的状态参数,,
Figure FDA0002616833170000043
表示网格单元j0中的状态参数;
Figure FDA0002616833170000044
表示网格单元i0相关的位置坐标,
Figure FDA0002616833170000045
表示网格单元j0相关的位置坐标;
针对i0和j0,利用所述量子随机存取存储器提取其对应的状态参数
Figure FDA0002616833170000046
Figure FDA0002616833170000047
并利用上述第二Oracle提取其对应的坐标信息
Figure FDA0002616833170000051
Figure FDA0002616833170000052
构建用于同时编码所述坐标信息
Figure FDA0002616833170000053
Figure FDA0002616833170000054
和所述状态参数
Figure FDA0002616833170000055
Figure FDA0002616833170000056
的第一子量子线路Oracle OA1,Oracle OA1用于实现:
Figure FDA0002616833170000057
构建编码
Figure FDA0002616833170000058
的第二子量子线路Oracle OA2,Oracle OA2用于实现:
Figure FDA0002616833170000059
构建第三子量子线路Oracle OA3,其中:所述第三子量子线路Oracle OA3是所述第一子量子线路Oracle OA3的转置共轭量子线路,Oracle OA3用于实现:
Figure FDA00026168331700000510
基于所述第一子量子线路Oracle OA1、所述第二子量子线路Oracle OA2和所述第三子量子线路Oracle OA3构建Oracle OA
8.根据权利要求5所述的方法,其特征在于,所述针对所有的网格单元,基于所述第三量子线路求解所述网格单元的线性系统方程,当所述网格单元的流体状态趋于稳定时,获得所述网格单元线性系统方程的状态参数表示的流体状态作为所述网格单元的目标状态,具体包括:
针对所有的网格单元,基于所述第三量子线路表示的所述线性系统方程,从所述网格单元的状态参数的初始值开始,对所述线性系统方程进行迭代直至bn趋近于0时,获取此时的ΔUn
根据获取的ΔUn对所述量子随机存取存储器中存储的Un进行更新,得到Un+1作为所述网格单元的目标状态。
9.根据权利要求8所述的方法,其特征在于,所述ΔUn为信息编码在量子态振幅上的模拟编码态数据;
所述根据ΔUn对所述量子随机存取存储器中存储的Un进行更新,得到Un+1作为所述网格单元的目标状态,具体包括:
构建第一数据转化量子线路Oracle Oconvert1,Oracle Oconvert1用于实现:
ΔUn|i>|0>→|i>|ΔUn>;
从所述量子随机存取存储器中提取Un
构建量子加法器线路Oracle Oadd,Oracle Oadd用于实现:
|i〉|Un〉|ΔU〉→|i〉|Un〉|Un+1>;
构建第二数据转化量子线路Oracle Oconvert2,Oracle Oconvert2用于实现:|i〉|Un>|Un+1>→|i>|0>|Un+1>;
获得由依次设置的所述第一数据转化量子线路Oracle Oconvert1、所述量子随机存取存储器、所述量子加法器线路Oracle Oadd和所述第二数据转化量子线路Oracle Oconvert2的量子线路编码的由量子态本征态编码的|i〉|0〉|Un+1〉作为所述线性系统方程所包含的状态参数目标值;
根据所述状态参数目标值得到所述网格单元的目标状态。
10.根据权利要求9所述的方法,其特征在于,所述根据所述状态参数目标值得到所述网格单元的目标状态,具体包括:
构建第三数据转化量子线路Oracle Oconvert3,Oracle Oconvert3用于实现:|i>|0>|Un+1>→Un+1|i〉;
利用所述第三数据转化量子线路Oracle Oconvert3将所述状态参数目标值转换为由量子态振幅编码的Un+1,获得各所述网格单元的目标状态。
11.一种基于量子算法的计算流体动力学模拟装置,其特征在于,所述装置包括:
网格量子线路构建模块,用于在利用有限体积法解析计算流体动力学的过程中,针对流体运动的离散化数值网格中各网格单元,构建表示所述网格单元坐标信息的第一量子线路、表示所述网格单元的状态参数的第二量子线路,其中,所述网格单元的状态参数存储在量子随机存取存储器中,所述量子随机存取存储器可操作处于量子叠加态的地址和数据;
参数量子线路构建模块,用于基于所述第一量子线路、所述第二量子线路、所述量子随机存取存储器构建表示所述网格单元流体状态变化的线性系统方程的参数的第三量子线路;
目标状态获取模块,用于针对所有的网格单元,基于所述第三量子线路求解所述网格单元的线性系统方程,当所述网格单元的流体状态趋于稳定时,获得所述网格单元线性系统方程的状态参数表示的流体状态作为所述网格单元的目标状态。
12.一种量子计算机设备,其特征在于,包括:量子线路、量子效应器及量子随机存取存储器,所述量子计算机设备在运行时实现权利要求1-10任一所述的基于量子算法的计算流体动力学模拟方法。
13.一种计算机可读存储介质,其特征在于,所述计算机可读存储介质内存储有计算机程序,所述计算机程序被处理器执行时实现权利要求1-10任一所述的基于量子算法的计算流体动力学模拟方法。
CN202010771592.2A 2020-08-04 2020-08-04 基于量子算法的计算流体动力学模拟方法、装置及设备 Active CN114091363B (zh)

Priority Applications (4)

Application Number Priority Date Filing Date Title
CN202010771592.2A CN114091363B (zh) 2020-08-04 2020-08-04 基于量子算法的计算流体动力学模拟方法、装置及设备
US18/017,571 US20230297745A1 (en) 2020-08-04 2020-12-29 Computational fluid dynamics simulation method and apparatus based on quantum algorithm, and device
PCT/CN2020/140833 WO2022027916A1 (zh) 2020-08-04 2020-12-29 基于量子算法的计算流体动力学模拟方法、装置及设备
EP20948742.0A EP4195090A4 (en) 2020-08-04 2020-12-29 METHOD AND DEVICE FOR COMPUTERIAL FLUID DYNAMICS SIMULATION BASED ON A QUANTUM ALGORITHM AND DEVICE

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010771592.2A CN114091363B (zh) 2020-08-04 2020-08-04 基于量子算法的计算流体动力学模拟方法、装置及设备

Publications (2)

Publication Number Publication Date
CN114091363A true CN114091363A (zh) 2022-02-25
CN114091363B CN114091363B (zh) 2023-08-08

Family

ID=80119960

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010771592.2A Active CN114091363B (zh) 2020-08-04 2020-08-04 基于量子算法的计算流体动力学模拟方法、装置及设备

Country Status (4)

Country Link
US (1) US20230297745A1 (zh)
EP (1) EP4195090A4 (zh)
CN (1) CN114091363B (zh)
WO (1) WO2022027916A1 (zh)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114745111A (zh) * 2022-04-12 2022-07-12 中南林业科技大学 基于键控链式受控非和硬币算子的量子短密钥密码方法
CN115809705A (zh) * 2022-11-21 2023-03-17 合肥本源量子计算科技有限责任公司 基于量子计算的流体动力学计算系统及量子计算机
WO2023231511A1 (zh) * 2022-05-30 2023-12-07 苏州元脑智能科技有限公司 一种量子数据加载方法、装置、设备和可读存储介质
WO2024007919A1 (zh) * 2022-07-06 2024-01-11 本源量子计算科技(合肥)股份有限公司 基于lbm的量子流动模拟方法、装置、介质及设备
CN118036764A (zh) * 2024-04-09 2024-05-14 国开启科量子技术(安徽)有限公司 用于量子相位编码的方法、装置、设备及存储介质

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113657014B (zh) * 2021-08-10 2024-08-23 南京师范大学 一种基于量子游走的临近空间大气状态模拟方法及装置
CN115130675B (zh) * 2022-09-02 2023-01-24 之江实验室 一种量子随机电路的多振幅模拟方法和装置
CN115905778B (zh) * 2022-11-15 2024-08-09 四川大学 一种基于量子计算的线性方程加速方法及系统
CN116306122B (zh) * 2023-03-07 2023-09-19 国家海洋环境预报中心 用于海洋数值模式的网格计算方法、计算机及存储介质
CN116228993B (zh) * 2023-05-08 2023-08-25 中国空气动力研究与发展中心计算空气动力研究所 一种网格边构建方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103414477A (zh) * 2013-07-12 2013-11-27 西安电子科技大学 一种构造量子卷积码状态转移图和网格图的方法
CN108154240A (zh) * 2017-12-29 2018-06-12 合肥本源量子计算科技有限责任公司 一种低复杂度的量子线路模拟系统
US20190197426A1 (en) * 2017-10-24 2019-06-27 Nippon Telegraph And Telephone Corporation Transformation apparatus, decision apparatus, quantum computation apparatus, and quantum machine learning system
CN111461335A (zh) * 2020-04-03 2020-07-28 合肥本源量子计算科技有限责任公司 基于mpi多进程的含噪声单量子逻辑门实现方法及装置

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103559740B (zh) * 2013-08-26 2016-06-29 空气动力学国家重点实验室 一种可实现交互操作的棱柱网格生成方法
CN104331271A (zh) * 2014-11-18 2015-02-04 李桦 用于cfd的并行计算方法及系统
US11493066B2 (en) * 2016-01-20 2022-11-08 Soliton Holdings Generalized jet-effect and enhanced devices
CN106780747B (zh) * 2016-11-30 2019-05-10 西北工业大学 一种快速分割cfd计算网格的方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103414477A (zh) * 2013-07-12 2013-11-27 西安电子科技大学 一种构造量子卷积码状态转移图和网格图的方法
US20190197426A1 (en) * 2017-10-24 2019-06-27 Nippon Telegraph And Telephone Corporation Transformation apparatus, decision apparatus, quantum computation apparatus, and quantum machine learning system
CN108154240A (zh) * 2017-12-29 2018-06-12 合肥本源量子计算科技有限责任公司 一种低复杂度的量子线路模拟系统
CN111461335A (zh) * 2020-04-03 2020-07-28 合肥本源量子计算科技有限责任公司 基于mpi多进程的含噪声单量子逻辑门实现方法及装置

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
A. FEDOSEYEV; T. BALD: ""Conductivity of quantum-dot crystal simulation using physics based models"", 《2011 NUMERICAL SIMULATION OF OPTOELECTRONIC DEVICES》 *
卢风顺; 陈波; 江雄: ""量子计算及其在空气动力学中的应用前景"", 《航空学报》 *
董果香: "" 基于量子流体动力学模型的半导体器件模拟"", 《电子设计工程》 *
魏世杰;王涛;阮东;龙桂鲁;: "量子算法的一些进展", 中国科学:信息科学, no. 10 *

Cited By (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114745111A (zh) * 2022-04-12 2022-07-12 中南林业科技大学 基于键控链式受控非和硬币算子的量子短密钥密码方法
CN114745111B (zh) * 2022-04-12 2024-04-30 中南林业科技大学 基于键控链式受控非和硬币算子的量子短密钥密码方法
WO2023231511A1 (zh) * 2022-05-30 2023-12-07 苏州元脑智能科技有限公司 一种量子数据加载方法、装置、设备和可读存储介质
WO2024007919A1 (zh) * 2022-07-06 2024-01-11 本源量子计算科技(合肥)股份有限公司 基于lbm的量子流动模拟方法、装置、介质及设备
CN115809705A (zh) * 2022-11-21 2023-03-17 合肥本源量子计算科技有限责任公司 基于量子计算的流体动力学计算系统及量子计算机
CN118036764A (zh) * 2024-04-09 2024-05-14 国开启科量子技术(安徽)有限公司 用于量子相位编码的方法、装置、设备及存储介质
CN118036764B (zh) * 2024-04-09 2024-08-06 国开启科量子技术(安徽)有限公司 用于量子相位编码的方法、装置、设备及存储介质

Also Published As

Publication number Publication date
CN114091363B (zh) 2023-08-08
EP4195090A1 (en) 2023-06-14
EP4195090A4 (en) 2024-06-26
US20230297745A1 (en) 2023-09-21
WO2022027916A1 (zh) 2022-02-10

Similar Documents

Publication Publication Date Title
CN114091363B (zh) 基于量子算法的计算流体动力学模拟方法、装置及设备
JP7471736B2 (ja) 量子系の基底状態エネルギーの推定方法、およびシステム
CN111599414B (zh) 一种基于量子计算机的全量子分子模拟方法
Cortese et al. Loading classical data into a quantum computer
Zulehner et al. Matrix-vector vs. matrix-matrix multiplication: Potential in DD-based simulation of quantum computations
CN114219076B (zh) 量子神经网络训练方法及装置、电子设备和介质
CN113273082A (zh) 具有异常块浮点的神经网络激活压缩
Penkovsky et al. Efficient design of hardware-enabled reservoir computing in FPGAs
CN114492814B (zh) 基于量子计算模拟目标体系能量的方法、装置及介质
WO2020221583A1 (en) System and method for molecular design on a quantum computer
Itani et al. Quantum algorithm for lattice Boltzmann (QALB) simulation of incompressible fluids with a nonlinear collision term
Mann et al. Development of a robust CNN model for capturing microstructure-property linkages and building property closures supporting material design
CN114492815B (zh) 一种基于量子化学计算目标体系能量的方法、装置及介质
CN114528996B (zh) 一种目标体系试验态初始参数的确定方法、装置及介质
Park et al. Convolution hierarchical deep-learning neural network (c-hidenn) with graphics processing unit (gpu) acceleration
CN115169565A (zh) 一种小分子化学体系的哈密顿量模拟方法和装置
Zhou et al. Quantum image scaling based on bilinear interpolation with decimals scaling ratio
CN114936646A (zh) 一种量子化数据处理方法和装置
CN115618663B (zh) 一种网格方程与物理方程耦合的量子求解方法及装置
Baker et al. Direct solution of multiple excitations in a matrix product state with block Lanczos
Pande A Comprehensive Review of Data Encoding Techniques for Quantum Machine Learning Problems
JP2002042104A (ja) 量子ソフトコンピューティングを使用した制御システムと制御方法
JP4758965B2 (ja) ユニタリ行列分解方法、ユニタリ行列分解装置、ユニタリ行列分解プログラム及び記録媒体
Hodapp et al. Equivariant Tensor Networks
Zhi et al. Learning odes via diffeomorphisms for fast and robust integration

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant