第十七講-并行計(jì)算

八、并行計(jì)算CIG-ComputationalInfrastructureforGeodynamics/GeoFEMProjecthttp://geofem.tokyo.rist.or.jp/QuakeSimProject/TeraShake–TypicalFramework/projects/scec/terashake/GEONGRIDProject:

常見的大型地球科學(xué)并行計(jì)算軟件系統(tǒng)GCNfileSourceCodeGeneratorPDEfileGIOfilePreprocessinginterfaceDatastructureMeshgenerationAllthesourcecodeforfiniteelementcomputinginterfacetosolversPostprocessingVisualizationtoolsinterfaceMeshgenerationtoolsComputingdomainconstructiontoolsMeshgenerationlibMeshpartitionlibDatapartitionlibSequentialandParallelsolverslibFinishedmodelslibAuxiliaryroutineslibAmira,GID,OpenDX,ParallelView,RocksclustersVizArchitectureprofileofparallelfiniteelementprogramgenerator/Computingenvironment張懷，2010pFEPG并行有限元程序自動(dòng)生成系統(tǒng)的組織架構(gòu)dispuvcoorxyfuncfunafunbfuncshap%1%2gaus%3mass%1load=fufv$c6pe=prmt(1)$c6pv=prmt(2)$c6fu=prmt(3)$c6fv=prmt(4)$c6fact=pe/(1.+pv)/(1.-2.*pv)funcfuna=+[u/x]funb=+[v/y]func=+[u/y]+[v/x]stifdist=+[funa;funa]*fact*(1.-pv)+[funa;funb]*fact*(pv)+[funb;funa]*fact*(pv)+[funb;funb]*fact*(1.-pv)+[func;func]*fact*(0.5-pv)*es,em,ef,Estifn,Estifv,*es(k,k),em(k),ef(k),Estifn(k,k),Estifv(kk),

goto(1,2),ityp1callseuq4g2(r,coef,prmt,es,em,ec,ef,ne)goto32callseugl2g2(r,coef,prmt,es,em,ec,ef,ne)goto33continueDOJ=1,NMATEPRMT(J)=EMATE((IMATE-1)*NMATE+J)EnddoPRMT(NMATE+1)=TIMEPRMT(NMATE+2)=DTprmt(nmate+3)=imateprmt(nmate+4)=numOtherelementmatrixcomputingSubsPDEexpressionContainsinformationofthephysicalmodel,suchasvariablesandequationsforgeneratingelementstiffnessmatrix.FortranSegmentscodesthatrealizethephysicalmodelatelementlevel.variablesequationAutomatedCodeGeneratorStep1:FromPDEexpressiontoFortransegmentsSegment1Segment2Segment3Segment4Step2:FromalgorithmexpressiontoFortransegments

doi=1,kdoj=1,kestifn(i,j)=0.0enddoenddodoi=1,kestifn(i,i)=estifn(i,i)doj=1,kestifn(i,j)=estifn(i,j)+es(i,j)enddoenddo

U(IDGF,NODI)=U(IDGF,NODI)*+ef(i)defistifSmassMloadFtypeemdtylstep0equationmatrix=[S]FORC=[F]SOLUTIONUwrite(s,unod)UendAlgorithmExpressionContainsinformationforformingglobalstiffnessmatrixforthemodel.FortranSegmentscodesthatrealizethephysicalmodelatgloballevel.StiffnessmatrixSegment5Segment6SUBROUTINEETSUB(KNODE,KDGOF,IT,KCOOR,KELEM,K,KK,*NUMEL,ITYP,NCOOR,NUM,TIME,DT,NODVAR,COOR,NODE,#SUBET.sub*U)implicitdoubleprecision(a-h,o-z)DIMENSIONNODVAR(KDGOF,KNODE),COOR(KCOOR,KNODE),*U(KDGOF,KNODE),EMATE(300),#SUBDIM.sub*R(500),PRMT(500),COEF(500),LM(500)#SUBFORT.sub#ELEM.subCWRITE(*,*)'ESEMEF='CWRITE(*,18)(EF(I),I=1,K)#MATRIX.subL=0M=0I=0DO700INOD=1,NNE………U(IDGF,NODI)=U(IDGF,NODI)#LVL.subDO500JNOD=1,NNE………500CONTINUE700CONTINUE………returnendProgramStencilFortranSegmentsgeneratedStep3:PlugFortransegmentsintoastencil,formingfinalFEprogramSegment1Segment2Segment4Segment3Segment5Segment6…………..Fromphysicalmodeltonumericalcomputation:

AdifferentapproachPartialDifferenceEquations(PDEs)PDE2FEMserverDataGrid(GEONandothers)

High-performancecomputingVianetworktouploadVianetworktouploadNumericalresultsVianetworktodownloadpost-processingAnalyzeandvisualitionAutomatedsourcecodegeneratorAutomaticsourcecodegeneratorfuncfuna=+[u/x]

………funf=+[u/y]+[v/x]

………dist=+[funa;funa]*d(1,1)+[funa;funb]*d(1,2)+[funa;func]*d(1,3)+[funb;funa]*d(2,1)+[funb;funb]*d(2,2)+[funb;func]*d(2,3)+[func;funa]*d(3,1)+[func;funb]*d(3,2)+[func;func]*d(3,3)+[fund;fund]*d(4,4)+[fune;fune]*d(5,5)+[funf;funf]*d(6,6)

load=+[u]*fu+[v]*fv+[w]*fw-[funa]*f(1)-[funb]*f(2)-[func]*f(3)-[fund]*f(4)-[fune]*f(5)-[funf]*f(6)PDEsCompletesourcecodeFEMModelingLanguageDataGrid(GEONandothers)PhysicalmodelModelresultsHPCCData=>???SWFSWF計(jì)算地球動(dòng)力學(xué)的學(xué)科交叉特色-機(jī)遇與挑戰(zhàn)并存并行有限元方法基礎(chǔ)理論研究及其自動(dòng)生成實(shí)現(xiàn)1997-2003：Lagrange乘子法，區(qū)域分解算法(DomainDecompositionMethod)，非連續(xù)變形分析算法(LagrangeMultiplierDiscontinuousDeformationAnalysis)的算法與并行自動(dòng)生成實(shí)現(xiàn);2002-2009：預(yù)條件Krylov子空間算法求解器的優(yōu)化改進(jìn)(Preconditionered

KrylovSubspaceIterativeSolvers)與并行有限元自動(dòng)生成技術(shù)的理論與算法實(shí)現(xiàn)；多場強(qiáng)耦合系統(tǒng)并行有限元求解算法；國內(nèi)空白，國際上目前只有類似的串行系統(tǒng)。國際上有類似有限元方法的自動(dòng)求解軟件包(Femlab,DiffPack,PDE2D等)，但是沒有自動(dòng)生成并行有限元源程序技術(shù)；已經(jīng)形成商業(yè)產(chǎn)品pFEPG(ParallelFiniteElementProgramGenerator)；計(jì)算地球動(dòng)力學(xué)若干模型的并行有限元計(jì)算算法研究與實(shí)現(xiàn)(2003-)板塊與活動(dòng)構(gòu)造（ActiveTectonics）-粘彈性并行有限元模型；區(qū)域性和全球地震波數(shù)值模擬-大規(guī)模并行顯式有限元算法與并行實(shí)現(xiàn);地幔對流與大陸演化過程-多場強(qiáng)耦合系統(tǒng)求解算法;基礎(chǔ)模型–線彈性模型與對流擴(kuò)散溫度場方程;計(jì)算地球動(dòng)力學(xué)基礎(chǔ)架構(gòu)中的若干輔助系統(tǒng)–跟蹤與科研需要(2005-)海量數(shù)據(jù)并行可視化設(shè)備(ParallelVisualization);

非結(jié)構(gòu)化網(wǎng)格與并行CVT/SCVT網(wǎng)格并行生成算法；Preliminaryresults通過高效的分區(qū)算法實(shí)現(xiàn)，用戶能自由指定分區(qū)，并自動(dòng)完成各個(gè)節(jié)點(diǎn)的負(fù)載均衡；1,200,000quadrilateralstructuredmeshpartition,sub-domainswithdifferentcolors(left)andsub-domainboundaries(right)50,000triangleunstructuredmeshpartition,(partitionedby10)研究探討一些地球科學(xué)問題中的基礎(chǔ)科學(xué)問題：1區(qū)域性活動(dòng)構(gòu)造-Maxwell粘彈性模型;2斷層系統(tǒng)的應(yīng)力演化與地震活動(dòng)性關(guān)系-粘-彈-塑性模型;3區(qū)域性板塊構(gòu)造演化歷史–強(qiáng)非線性的非牛頓流體熱演化問題、大變形（幾何非線性）；4地球自轉(zhuǎn)快速變化后的全球應(yīng)力場變化分析；5區(qū)域性應(yīng)力演化歷史與地震活動(dòng)性分析；關(guān)系到國際民生問題的研究-大的工程應(yīng)用方面的基礎(chǔ)研究工作6二維、三維地震波強(qiáng)地面運(yùn)動(dòng)數(shù)值模擬；7海嘯數(shù)值模擬，危險(xiǎn)性分析；8港口安全性評估；基礎(chǔ)設(shè)置建設(shè)其它94D海量數(shù)據(jù)體的并行可視化和體繪制的并行顯示；

10GPU計(jì)算與GPU+CPU層次結(jié)構(gòu)計(jì)算的算法研究；

10建立計(jì)算地球動(dòng)力學(xué)網(wǎng)格計(jì)算基礎(chǔ)架構(gòu)（國際合作，跟蹤研究）；地球動(dòng)力學(xué)方面的科學(xué)研究內(nèi)容

能不能做一些關(guān)系國計(jì)民生的問題？地震海嘯港口安全滑坡泥石流颶風(fēng)……地震災(zāi)害預(yù)測與救助指導(dǎo)E-Science系統(tǒng)現(xiàn)有工作基礎(chǔ)：全球地震波傳播數(shù)值模擬；區(qū)域地震波數(shù)值模擬-采用有限元方法，具有差分法和譜元法不具備的優(yōu)勢，同時(shí)采用合理的并行顯格式算法，同樣可以達(dá)到差分法的效率和譜元法的高分辨率；突出的優(yōu)點(diǎn)是可以考慮各種類型的單元；PeakGroundAcceleration(PGA)分析；功率譜密度（PowerSpectralDensity）分析，確定區(qū)域本征值分布特征；非結(jié)構(gòu)化網(wǎng)格生成技術(shù)-可以充分考慮多尺度問題；數(shù)據(jù)庫累積；并行可視化系統(tǒng)；網(wǎng)格國際合作，整合即可完成本研究項(xiàng)目；PFEPG并行有限元程序自動(dòng)生成系統(tǒng)并行程序與串行程序的區(qū)別在于：1.是否分區(qū)，2.是否存在主從節(jié)點(diǎn)，3.是否存在節(jié)點(diǎn)之間的數(shù)據(jù)交換。PFEPG工作原理圖客戶端——PFEPG服務(wù)器——并行機(jī)三者之間的關(guān)系應(yīng)用pFEPG進(jìn)行并行計(jì)算：

客戶機(jī)不是并行機(jī)的節(jié)點(diǎn)pFEPG并行程序結(jié)構(gòu)解析算法基礎(chǔ)：區(qū)域分解求解問題范圍：動(dòng)態(tài)問題穩(wěn)態(tài)問題線性問題非線性問題元件程序組成（六個(gè)模塊）前處理分區(qū)模塊（partition）START、BFT、Solve、E和U元件程序注：其中E和U元件程序根據(jù)用戶給出的偏微分方程表達(dá)式和算法表達(dá)式由系統(tǒng)自動(dòng)生成，其余幾個(gè)元件程序由系統(tǒng)給出，不隨表達(dá)式的變動(dòng)而變動(dòng)。本系統(tǒng)的基本思想是元件化程序設(shè)計(jì)和程序自動(dòng)生成技術(shù)相結(jié)合，根據(jù)區(qū)域分解法的特點(diǎn)，將有限元程序分成：前處理分區(qū)、START、BFT、Solve、E和U等六個(gè)模塊，其中E和U元件程序程序根據(jù)用戶給出的偏微分方程表達(dá)式和算法表達(dá)式由系統(tǒng)自動(dòng)生成，其余幾個(gè)元件程序由系統(tǒng)給出，不隨表達(dá)式的變動(dòng)而變動(dòng)。PFEPG并行有限元程序自動(dòng)生成系統(tǒng)前處理分區(qū)模塊：該模塊又分成四個(gè)元件程序：partition、getpart、spart、sgetpart元件；該模塊的功能事根據(jù)用戶提供的前處理數(shù)據(jù)將整個(gè)求解區(qū)域根據(jù)用戶的要求合理劃分成數(shù)個(gè)子區(qū)域。前處理分區(qū)模塊主進(jìn)程：Partition.fGetpart.fPartition.dat從進(jìn)程：Spart.fSgetpart.f該功能模塊是根據(jù)用戶提供的前處理數(shù)據(jù)將整個(gè)求解區(qū)域根據(jù)用戶的要求合理劃分成數(shù)個(gè)子區(qū)域，partition給出各個(gè)子區(qū)域所包含的整體節(jié)點(diǎn)號，局部節(jié)點(diǎn)號和整體節(jié)點(diǎn)編號的對應(yīng)關(guān)系，并將各子區(qū)域的一些前處理信息，包括節(jié)點(diǎn)數(shù)、時(shí)間信息等，以及各區(qū)域所含節(jié)點(diǎn)的坐標(biāo)值分塊發(fā)送給與各子區(qū)域?qū)?yīng)的從進(jìn)程。Getpart根據(jù)partition提供的各子區(qū)域所含節(jié)點(diǎn)信息得到各自的局部單元信息、初邊值信息、節(jié)點(diǎn)ID信息以及各子區(qū)域節(jié)點(diǎn)的各自由度對應(yīng)的整體方程號等，并發(fā)送各相應(yīng)從進(jìn)程。Spart元件接收來自主進(jìn)程partition元件的坐標(biāo)值；Sgetpart接收來自主進(jìn)程getpart元件發(fā)送的消息。前處理分區(qū)模塊有限元網(wǎng)格分區(qū)示例分區(qū)圖START元件：功能：在各子區(qū)域內(nèi)給出每個(gè)節(jié)點(diǎn)的各個(gè)自由度與將要形成的代數(shù)方程組的變量（方程號）的對應(yīng)關(guān)系。得到整體剛度矩陣的存儲格式，該元件程序只在從進(jìn)程執(zhí)行。BFT元件：功能：計(jì)算邊值，更新時(shí)間。該元件由兩個(gè)程序：bftm和bfts組成，分別由主進(jìn)程和從進(jìn)程執(zhí)行，bftm由主進(jìn)程執(zhí)行，主要是給出每一時(shí)刻解的邊值，對時(shí)間進(jìn)行更新，進(jìn)行結(jié)束與否的判斷，并將判斷結(jié)果和新的時(shí)間和各區(qū)域解的邊值發(fā)送給從進(jìn)程

；bftm.fbfts由從進(jìn)程執(zhí)行，接收由主進(jìn)程發(fā)送的新的時(shí)間和解的邊值。bfts.fE元件：功能：該元件程序用于單元計(jì)算，包括單剛，單質(zhì)和單元載荷計(jì)算，并把它們由節(jié)點(diǎn)自由度表示轉(zhuǎn)換成代數(shù)方程組的變量表示，同時(shí)處理邊界約束條件，并形成代數(shù)方程組的系數(shù)矩陣和右端項(xiàng)。SOLV元件：功能：組裝和求解整體線性代數(shù)方程組。由兩個(gè)程序組成：solv和msolvsolv由從進(jìn)程執(zhí)行Msolv由主進(jìn)程執(zhí)行U元件：功能：方程組解轉(zhuǎn)換回節(jié)點(diǎn)自由度值，并進(jìn)行其它后處理計(jì)算包括u元件程序和mu元件程序U元件由從進(jìn)程執(zhí)行，進(jìn)行解的轉(zhuǎn)換和后處理，并將結(jié)果傳輸?shù)街鬟M(jìn)程Mu元件程序由主進(jìn)程執(zhí)行，接收來自于各子區(qū)域的解，整理得到整體區(qū)域上的解，更新子區(qū)域公共邊界解，結(jié)果存盤處理。為了使并行程序運(yùn)行結(jié)果與串行程序運(yùn)行結(jié)果保持一致，同時(shí)避免由于并行I/O所帶來的程序設(shè)計(jì)上的復(fù)雜性，我們規(guī)定從進(jìn)程不參與磁盤文件的I/O，所有磁盤文件的I/O操作均由主進(jìn)程完成。需要保存的數(shù)據(jù)寫入預(yù)先開辟的內(nèi)存段中，實(shí)現(xiàn)數(shù)據(jù)文件與特定內(nèi)存段之間的對應(yīng)關(guān)系。這些數(shù)據(jù)文件保存都將是一些數(shù)組，單個(gè)的數(shù)據(jù)被作為公共變量在各元件程序之間進(jìn)行傳遞。數(shù)據(jù)結(jié)構(gòu)眾所周知，并行程序的調(diào)試非常困難，PFEPG系統(tǒng)的這種工作模式最大程度上降低了程序的調(diào)試難度。用戶上傳的有限元表達(dá)式和由服務(wù)器生成的計(jì)算源程序均是文本文件的形式，其大小一百K以內(nèi)，在Internet上傳輸速度很快，一般來說只需十幾秒即可生成用戶所需的全部源程序，并傳輸?shù)接脩舻谋镜貦C(jī)。下面說明PFEPG系統(tǒng)的使用步驟：1．填寫Partition文件Partition.dat文件中是一些前處理分區(qū)程序需要的常數(shù)，內(nèi)容一般如下41221其中，第一行的第一個(gè)整數(shù)4表示總的分塊數(shù)，第二個(gè)整數(shù)1表示采用的分區(qū)算法（1表示采用圖論的分區(qū)算法，2表示利用坐標(biāo)分區(qū)的算法）；第二行的第一個(gè)2表示x方向剖分?jǐn)?shù)，第二個(gè)2表示y方向剖分?jǐn)?shù)，1表示z方向剖分?jǐn)?shù)。第一行的第一個(gè)數(shù)是第二行三個(gè)數(shù)的乘積，在做大規(guī)模計(jì)算時(shí)，需要修改partition.dat文件來增加分塊數(shù)，但這個(gè)關(guān)系需要保持。如果用戶不填寫此文件，系統(tǒng)會(huì)產(chǎn)生如上的缺省文件。2．用戶填寫完公式文件（PDE,GCN,GIO等）并通過界面?zhèn)魉偷椒?wù)器后，運(yùn)行mpi命令如下：mpigcn文件的前綴名即可產(chǎn)生并行有限元源程序，并傳回客戶機(jī)。數(shù)據(jù)文件說明partition.dat功能：指定和記錄計(jì)算時(shí)分區(qū)數(shù)目和方法格式：41221M：總的分區(qū)數(shù)(>=2)N：分區(qū)方法，只可取1,2，分別表示按圖論的分塊方法和按坐標(biāo)的分塊方法MxMyMz:三個(gè)坐標(biāo)軸方向的分區(qū)數(shù)M=Mx*My*Mz…MN…MxMyMz3．把以下文件傳送到并行機(jī)全部并行fortran源程序（*.f文件），fegen.b、runmpi、time0、nbefile、partition.dat、*.io。如果使用的客戶端是Window系統(tǒng)下的，對于上傳的文件，一般需要在并行機(jī)上運(yùn)行如下命令dos2unix*.*將上傳的文件轉(zhuǎn)成Unix文件格式。對于由多臺PC機(jī)或工作站組成的并行機(jī)，PFEPG客戶端安裝在其中一臺PC或工作站上的情況，不需要傳送任何文件。4．在并行機(jī)上運(yùn)行如下命令csh–xfegen.b通過編譯器產(chǎn)生并行的可執(zhí)行程序。5．運(yùn)行前處理，mpigid，cshpreproc.mpi。6．在并行機(jī)上運(yùn)行如下命令csh–xrunmpi進(jìn)行并行計(jì)算。7．運(yùn)行后處理，cshpostproc.mpi。pFEPG程序代碼淺析pFEPG封裝的幾個(gè)基本通信函數(shù)sendint(des,src,idata)：從源進(jìn)程src發(fā)送一個(gè)整型數(shù)idata到目的進(jìn)程dessendr(des,src,rdata)：從源進(jìn)程src發(fā)送一個(gè)實(shí)型數(shù)rdata到目的進(jìn)程dessendai(des,src,idata,len)：從源進(jìn)程src發(fā)送一個(gè)長度為len的整型數(shù)組idata到目的進(jìn)程dessendar(des,src,rdata,len)：從源進(jìn)程src發(fā)送一個(gè)長度為len的實(shí)型數(shù)組rdata到目的進(jìn)程despFEPG程序代碼淺析pFEPG封裝的幾個(gè)基本通信函數(shù)recvint(des,src,idata)：從源進(jìn)程src接收一個(gè)整型數(shù)idata到目的進(jìn)程desrecvr(des,src,rdata)：從源進(jìn)程src接收一個(gè)實(shí)型數(shù)rdata到目的進(jìn)程desrecvai(des,src,idata,len)：從源進(jìn)程src接收一個(gè)長度為len的整型數(shù)組idata到目的進(jìn)程desrecvar(des,src,rdata,len)：從源進(jìn)程src接收一個(gè)長度為len的整型數(shù)組rdata到目的進(jìn)程despFEPG程序代碼淺析代碼中的幾個(gè)固定變量和數(shù)組說明幾個(gè)常用的固定變量：numblk：分塊數(shù)/從進(jìn)程個(gè)數(shù)nsource：=0，主進(jìn)程進(jìn)程號iblk：塊號/從進(jìn)程進(jìn)程號主要用于消息通信pFEPG程序代碼淺析代碼中的幾個(gè)固定變量和數(shù)組說明幾個(gè)常用的固定數(shù)組：knode_iblk(numblk)：子區(qū)域節(jié)點(diǎn)總數(shù)n_update(numblk)：子區(qū)域內(nèi)部節(jié)點(diǎn)數(shù)mglnod(knode)：整體節(jié)點(diǎn)編號對應(yīng)的子區(qū)域節(jié)點(diǎn)編號mlgnod(numblk,max_lknode)：子區(qū)域節(jié)點(diǎn)編號對應(yīng)的整體節(jié)點(diǎn)編號iorder(numblk,max_lknode)：子區(qū)域節(jié)點(diǎn)標(biāo)識，1為內(nèi)部節(jié)點(diǎn)，0為擴(kuò)展節(jié)點(diǎn)主要用于解的轉(zhuǎn)換和整理算例：三維熱固耦合問題

熱固耦合問題模型控制方程熱傳導(dǎo)方程彈性體方程其中為彈性體位移，為作用于彈性體上的體積力。分別表示應(yīng)力和應(yīng)變向量。A是一個(gè)矩陣微分算子D為彈性矩陣Q為溫度場產(chǎn)生的熱應(yīng)力向量，由下式確定：原型文件et.gcnet.gioheat.vdesdisp.vdessdisp.vde客戶端需要準(zhǔn)備的文件heat.vdedispucoorxyzshap%1%2gaus%3mateekq1.0d01;0.0d0;vectxxyzstifdist=+[u/x_i;u/x_i]*ekload=+[u]*qendsdisp.vdeFUNC$c6fact=pe/(1.+pv)/(1.-2.*pv)$c6ft=(1.+pv)*fact*alfaex=+[u/x]

ey=+[v/y]

ez=+[w/z]eyz=+[w/y]+[v/z]exz=+[w/x]+[u/z]exy=+[u/y]+[v/x]

stifdist=+[em_i;em_j]*de_i_j*fact+[en_i;en_i]*fact*(0.5-pv)

load=+[u_i]*fd_i+[em_i]*fe_i

enddispu,v,wcoorx,y,zfuncexeyezeyzexzexyshap%1%2gaus%3coeftnmatepepvalfafxfyfz1.0e9;0.3;1.0e-5;0;0;7e3;vectemexeyezvecteneyzexzexyvectuuvwvectfdfxfyfzvectfeft*tnft*tnft*tnmatrde33(1.-pv)pvpvpv(1.-pv)pvpvpv(1.-pv)ssdisp.vdedispsx,sy,sz,syz,sxz,sxycoorx,y,zcoefu,v,w,tnshap%1%2gaus%3mass%11.0matepepvalfafxfyfz1.0e9;0.3;1.0e-5;0;0;7e3;vectfssxsyszsyzsxzsxyvectfffsxfsyfszfsyzfsxzfsxyvectemexeyezvecteneyzexzexyvectfmfsxfsyfszvectfnfsyzfsxzfsxyvectfefttfttfttmatrde33(1.-pv)pvpvpv(1.-pv)pvpvpv(1.-pv)stif$c6fact=pe/(1.+pv)/(1.-2.*pv)$cvftt=(1.+pv)*alfa*tn*fact$cvex=+{u/x}$cvey=+{v/y}$cvez=+{w/z}$cveyz=+{w/y}+{v/z}$cvexz=+{w/x}+{u/z}$cvexy=+{u/y}+{v/x}$cvfm_i=+fact*de_i_j