基于分数阶微积分的Kelvin_Voigt流变模型_郭佳奇
求解三维空间分数阶对流扩散方程的Douglas-Gunn格式

!第"#卷第#期郑州大学学报!理学版"$%&’"#(%’#!)*#+年,月-./012340%56278.!(9:.;<7.=>."?9@.)*#+收稿日期!)*#AB*#B *C 基金项目!国家自然科学基金项目!##ED#)C)".作者简介!聂玉峰!#+CA %"#男#陕西西安人#教授#主要从事高性能数值计算方法’计算材料学’并行计算研究#=B I 97&&M \271N 2[L5.1>5.<2$通信作者&胡嘉卉!#+A*%"#女#河南郑州人#博士研究生#主要从事偏微分方程数值解研究#=B I 97&&05R 0NI 97&.2[L5.1>5.<2.求解三维空间分数阶对流扩散方程的!9D )?4>=U D ::格式聂玉峰#!!胡嘉卉#!)!!王俊刚#!#.西北工业大学计算科学研究中心!陕西西安D#*#)+$).河南工业大学理学院!河南郑州E"***#"摘要!由于分数阶导数的非局部性特征#在模拟反常扩散现象时使用分数阶偏微分方程具有更好的效果#但是分数阶导数的非局部性也给数值分析和计算带来了很大困难#尤其在多维空间情形下.通过对经典d%53&9S B V 522格式的推广#提出一种求解三维空间分数阶对流扩散方程!SL9<1\@9<:7%29&9>81<:7%2>7\\5S 7%21O59:7%2#;F G d ="的交替方向隐!9&:1@29:723>7@1<:7%27I L&7<7:#G d ^"差分格式#并用矩阵法证明了其稳定性和收敛性.用数值算例进一步验证了该格式在空间和时间方向均具有较高的二阶收敛精度#可以高效地求解三维;FG d =.关键词!三维;F G d =$G d ^格式$Y @92JB (7<%&S %2格式$d %53&9S B V 522格式$稳定性$收敛性中图分类号!e)E#文献标志码!G 文章编号!#CD#B CAE#!)*#+"*#B **EEB *D !"#!#*’#,D*"Q R .7S S 2.#CD#B CAE#’)*#A*##$%引言近年来#自然界中的反常扩散现象受到科研人员的广泛关注#为研究其独特的物理过程#常常用分数阶偏微分方程建立相应的数学模型.其中#在包含对流和超扩散两个物理过程的散布现象中#粒子束的传播与经典的布朗运动模型不再一致#此时把经典对流扩散方程中的空间二阶导数替换成分数阶导数构建的空间分数阶对流扩散方程!;FG d ="能更准确地模拟这一现象.分数阶导数或积分具有非局部性#这给相应方程的求解带来了很大困难#在大多数情况下很难得到解析解#因此研究可靠而有效的数值方法就显得尤为重要.目前常用的数值解法包括有限差分法(#X ))’有限元法(,X E )’有限体积法(")’配点法(C )以及谱方法(D )等.由于三维模型在科学研究中有广泛应用#本文考虑有限区域上带有零d 7@7<0&1:边界条件的三维;F G d =的数值求解.分数阶导数是一个非局部算子#这就使得离散;FG d =得到的线性系统的刚度矩阵不再是稀疏矩阵#导致计算工作量和存储量都非常大#尤其在多维空间情形下.目前求解三维;FG d =的数值方法还比较少.文献(A )采用了一种交替方向稳定法!9&:1@29:723>7@1<:7%27I L&7<7:I 1:0%>#G d ^"差分格式求解三维;F G d =#并提高了精度.文献(+)提出了一种求解三维分数阶扩散方程的G d ^差分格式.文献(#*)研究了一种求解三维空间分数阶扩散方程的快速迭代G d ^有限差分方法.在本文中#我们将提出一种求解三维;F G d =的有效的G d ^有限差分格式#这种方法是将经典的d %53&9S B V 522格式中的二阶中心差分算子推广为包含左’右分数阶导数离散算子及一阶中心差分算子在内的复杂算子得到的#同时给出了该格式的稳定性和收敛性的必要证明.最后用数值实验验证了理论分析的结果.&%三维F ’7!L 及其!9D )?4>=U D ::格式本文考虑三维;FG d =及它的初边值条件为-@!I #C #J #""-"-:#I Q R "I @!I #C #J #""#:)I R "I +@!I #C #J #""#P #C Q R !C@!I #C #J #""# Copyright©博看网 . All Rights Reserved.!第#期聂玉峰#等$求解三维空间分数阶对流扩散方程的d%53&9S B V 522格式P )C R !C +@!I #C #J #""#2#J QR %J @!I #C #J #""#2)J R %J +@!I #C #J #""#X #-@!I#C #J #""-I#X )-@!I #C #J #""-C #X ,-@!I#C #J #""-J#!!I #C #J #""#!I #C #J #""(*/!*#Z )#!#"@!I #C #J #""-*#!!I #C #J #""(-*/(*#Z )#!)"@!I #C #J #*"-@*!I #C #J "#!!I #C #J "(*#!,"其中&#K "#!#%K )$*-!I Q #I +"/!C Q #C +"/!J Q #J +"4S ,$:##:)#P ##P )#2##2)!*是,个空间方向的左’右扩散系数$X #’X )’X ,分别是,个空间方向的对流系数$!!I #C #J #""是源项.方程!#"中的分数阶导数是_71I 922B ]7%587&&1型的#即函数A!I "的"!#K "K )"阶_71I 922B ]7%587&&1左导数和右导数分别定义为I Q R "IA !I "-#’!).""-)-I ),I I Q!I.4"#."A !4">4和I R "I +A !I "-#’!).""-)-I ),I +I!4.I "#."A !4">41&,&%S (<@4::=V (9D Q(??<分数阶导数的离散设*#’*)’*,和,为正整数#3#-!I +.I Q "L *##3)-!C +.C Q "L *)#3,-!J +.J Q "L *,#.-Z L ,分别是一致的空间步长和时间步长#由此定义的空间和时间的剖分为I ;-I Q #;3#;-*###-#*##C G -C Q #G 3#G -*###-#*)#J 9-J Q #93#9-*###-#*,#"7-7.#7-*###-#,1设@7;#G #9表示@!I ;#C G #J 9#"7"的近似值#!7;#G #9-!!I ;#C G #J 9#"7"1采用文献(##)中的方法离散方程!#"中的分数阶导数1以I方向为例#有IQR "I@!I #C G #J 9#"7"I -I ;-#’!E .""3"##;##5-*?"5@!I ;.5###C G #J 9#"7"#_!3)#"#IR "I +@!I #C G #J 9#"7"I -I ;-#’!E .""3"##*#.;##5-*?"5@!I ;#5.##C G #J 9#"7"#_!3)#"#其中系数为"5-##5-*#.E #),."#5-##C .)"."#,,."#5-)#!5##",.".E 5,."#C !5.#",.".E !5.)",."#!5.,",."#5!,1记)#"#I @7;#G #9-#’!E .""3"##;##5-*?"5@7;.5###G #9#!E ")."#I @7;#G #9-#’!E .""3"##*#.;##5-*?"5@7;#5.##G #91!""同时#用中心差分近似对流项的一阶导数#记R "#I @7;#G #9-!@7;###G #9.@7;.##G #9"L )3#1!C "为了便于表示#进一步引入记号R ^"#I @7;#G #9-X #R "#I @7;#G #9#)^#"#I @7;#G #9-:#)#"#I @7;#G #9#)^."#I @7;#G #9-:))."#I @7;#G #91C 和J 方向的记号可以类似表示.&,+%三维F ’7!L 的有限差分近似用公式!E "i !C "离散空间导数#时间方向采用Y@92JB (7<%&S %2格式.记)"#I o -)^#"#I #)^."#I #R ^"#I #)!#C o -)^#!#C #)^.!#C #R ^!#C #)%#J o -)^#%#J #)^.%#J #R ^%#J #然后方程!#"就可以表示为!#..))"#I ..))!#C ..))%#J "@!I ;#C G #J 9#"7##"-!##.))"#I #.))!#C #.))%#J "@!I ;#C G #J 9#"7"#.!!I ;#C G #J 9#"7##L )"#+7##;#G #91!D "存在正常数$^#使得+7##;#G #9’$^.!.)#3)##3))#3),"1接下来#在方程!D "中用近似值@7;#G #9代替函数值@!I ;#C G #J 9#"7"#并去掉高阶项#得到方程!#"的全离散格式#"E Copyright©博看网 . All Rights Reserved.郑州大学学报!理学版"第"#卷!#..))"#I ..))!#C ..))%#J "@7##;#G #9-!##.))"#I #.))!#C #.))%#J "@7;#G #9#.!7##L );#G #91!A "!!为便于计算#下面构造G d ^差分格式.将方程!A "左边加上高阶项!.)E )"#I )!#C #.)E )"#I )%#J #.)E )!#C )%#J ",!@7##;#G #9.@7;#G #9"..,A)"#I )!#C )%#J !@7##;#G #9#@7;#G #9"#再把适当的部分移项到方程的右边并分解因式#得!#..))"#I "!#..))!#C "!#..))%#J "@7##;#G #9-!##.))"#I "!##.))!#C "!##.))%#J "@7;#G #9#.!7##L );#G #91!+"采用经典的d %53&9S B V 522格式分解式!+"得到G d ^格式#即!#..))"#I "!@+-!.)"#I #.)!#C #.)%#J "@7;#G #9#.!7##L );#G #9$!#*"!#..))!#C"!@++-!@+$!##"!#..))%#J "!@-!@++$!#)"!@-@7##;#G #9.@7;#G #91!#,"与经典的d %53&9S B V 522格式不同#,个方向的二阶中心差分算子在此处分别被替换为)"#I ’)!#C ’)%#J #它们是包含了左’右分数阶导数离散算子等在内的复杂算子#可以认为是经典d %53&9S B V 522格式在求解分数阶方程中的推广.接下来#我们将给出收敛性和稳定性的必要证明.+%收敛性和稳定性分析显然#如果在格式!#*"i !#,"中消去中间解变量#则得到格式!+"#即格式!#*"i !#,"和!+"是等价的.下面用矩阵法证明格式!+"是无条件稳定和收敛的.首先把方程!+"表示成矩阵形式#令"7-(@7######@7)#####-#@7*#.######@7##)###@7)#)###-#@7*#.##)###-#@7##*).####@7)#*).####-#@7*#.##*).####@7####)#@7)###)#-#@7*#.####)#-#@7##*).##)#@7)#*).##)#-#@7*#.##*).##)#-#@7##*).##*,.##@7)#*).##*,.##-#@7*#.##*).##*,.#)P #!#E "为了书写简单起见#我们用记号"7-!!(@7;#G #9);-##)#-#*#.#"G -##)#-#*).#"9-##)#-#*,.#表示式!#E "#类似的表示还有87-!!(!7;#G #9);-##)#-#*#.#"G -##)#-#*).#"9-##)#-#*,.#1记0^I -:#.)’!E .""3"#71710"#:).)’!E .""3"#71710P"#X #.E 3#71711#!#""0^C -P #.)’!E .!"3!)710!17#P ).)’!E .!"3!)710P!17#X ).E 3)71117#!#C "0^J -2#.)’!E .%"3%,0%1717#2).)’!E .%"3%,0P%1717#X ,.E 3,11717#!#D "其中0"#0!#0%是P %1L&7:4矩阵#0"和1分别表示为0"-"#?"**-**")?"#?"**-*",?")?"#?"*-*222222"*#.)---?"#?"*"*#.#?"*#.)?"*#.,-?")?"##1-*#*-**.#*#*-**.#*#-*222222***-*#***-.#*#其中&7是单位矩阵$符号1表示b @%21<J1@积(#)).0!’0%与0"类似.利用上述记号#式!+"可以写为!7.0^I "!7.0^C "!7.0^J""7##-!7#0^I "!7#0^C "!7#0^J""7#.87##L )1!#A "!!为了证明式!#A "的稳定性和收敛性#下面列出一些相关的引理和定理.CE Copyright©博看网 . All Rights Reserved.!第#期聂玉峰#等$求解三维空间分数阶对流扩散方程的d%53&9S B V 522格式引理&(#,)!一个7阶实矩阵0是正定的#当且仅当矩阵,-!0#0P"L )是正定的$,是正定的#当且仅当,的特征值都是正的.引理+(#,)!设0是一个7阶复矩阵#0T 表示H 的共轭转置#记,-!0#0T "L)#则对于0的任意特征值(#它的实部满足不等式(I 72!,"’+!(!0""’(I 9‘!,"#这里(I 72!,"和(I 9‘!,"分别表示,的最小和最大特征值.定理&(A )!设0"是式!#""中的P %1L&7:4矩阵#则对于0"的任意特征值(#有+!(!0"""K *#并且0"是负定矩阵1同时#+!(!=#0"#=)0P """K *#=##=)!*#=)##=))8*1引理-!设0(S 9/7#1(S 8/B #2(S T /[#则!011"12-01!112"1证明!此结论可以由b @%21<J1@积的定义直接得到.引理.!设0#1(S 9/7#2(S B /"#则有!0#1"12-012#112#21!0#1"-210#2111证明!此结论可以由b @%21<J1@积的定义直接得到.引理/(#))!设0(S 9/7#1(S 8/B #2(S 7/T #3(S B /"#则!011"!213"-02113!(S 98/T ""1引理0(#))!对于任意的矩阵0和1#有!011"P -0P 11P .引理1(#))!设矩阵0(S 7/7有特征值0(;17;-##矩阵1(S 9/9有特征值01G 19G -#1则矩阵011的97个特征值为(#1##-#(#19#()1##-#()19#-#(71##-#(7191为了叙述并证明下述引理和定理#记%,%表示矩阵的)B 范数.引理2%设0(S 7/7是正定矩阵#则对任意的+(S 且+]*#有%!7#+0".#%K #1证明!由矩阵)B 范数的定义#有%!7#+0".#%)-I 9‘I 8*!!7#+0".#I #!7#+0".#I "!I#I "1设C -!7#+0".#I #则有%!7#+0".#%)-I9‘C 8*!C #C "!!7#+0"C #!7#+0"C "-#I 72C 8*(##)+!0C #C "!C #C "#+)!0C #0C "!C#C ")K #1引理6%设0(S 7/7是正定矩阵#则对任意的+(S 且+]*#有%!7#+0".#!7.+0"%K #1证明!由矩阵)B 范数的定义并记C-!7#+0".#I #可得%!7.+0"!7#+0".#%)-I 9‘I 8*!!7.+0"!7#+0".#I #!7.+0"!7#+0".#I "!I #I "-I 9‘C 8*!!7.+0"C #!7.+0"C "!!7#+0"C #!7#+0"C "-I 9‘C 8*!C #C ".)+!0C #C "#+)!0C #0C "!C #C "#)+!0C #C "#+)!0C #0C "K #1!!引理&$%设0#1#7(S 7/7#0和1乘积可交换#且!7.0".#’!7.1".#存在#则!7#0"与!7.1".##!7.0".#与!7.1".#也是乘积可交换的.证明!首先#由01-10#不难验证!7a 0"!7.1"-!7.1"!7a 0"1所以有!7#0"!7.1".#-!7.1".#!7.1"!7#0"!7.1".#-!7.1".#!7#0"!7.1"!7.1".#-!7.1".#!7#0"#!7.0".#!7.1".#-!!7.1"!7.0"".#-!!7.0"!7.1"".#-!7.1".#!7.0".#1定理+%由式!#""i !#D "定义的矩阵0^I ’0^C ’0^J 是负定的.证明!记0I -:#.)’!E .""3"#0"#:).)’!E .""3"#0P"#X #.E 3#1#则有0PI-:#.)’!E .""3"#0P"#:).)’!E .""3"#0"#X #.E 3#1P #0I #0P I )-!:##:)".)’!E .""3"#0"#0P")#由定理##可得+!(!0I #0P I )""K *#也就是(!0I #0P I )"K *#而且由引理,#E 和C #有0^I #0^PI)-7171!0I #0PI)"#再根据引理D #可以得到(!0^I #0^PI)"K *1最后由引理)和引理#可得+!(!0^I""K *#且0^I 是负定的1类似地#可以证明0^C 和0^J 是负定的.定理-%差分格式!+"是无条件稳定的.DE Copyright©博看网 . All Rights Reserved.郑州大学学报!理学版"第"#卷证明!我们利用差分格式!+"的矩阵形式!#A "证明.设"7和#7分别是格式!#A "对应于初值"*和#*的解#记’7o -"7.#7-!!(’7;#G #9);-##)#-#*#.#"G -##)#-#*).#"9-##)#-#*,.#1由式!#A "可得’7满足方程!7.0^I "!7.0^C "!7.0^J"’7##-!7#0^I "!7#0^C "!7#0^J"’71!#+"!!矩阵0^I ’0^C ’0^J 是乘积可交换的.事实上#只需验证!#""i !#D "中b @%21<J1@积形式的矩阵都是乘积可交换的.由引理,和"#可得!71710""!710!17"-!71!710"""!71!0!17""-71!!710""!0!17""-71!0!10""#同时!710!17"!71710""-!71!0!17""!71!710"""-71!!0!17"!710"""-71!0!10""#等等1由0^I ’0^C ’0^J 乘积可交换并利用引理#*#方程!#+"可以写为’7-!!7.0^J".#!7#0^J ""7!!7.0^C".#!7#0^C ""7!!7.0^I".#!7#0^I""7’*1由定理)和引理+可知%!7.0^I ".#!7#0^I "%K ##所以!7.0^I ".#!7#0^I "的谱半径小于##因此当7)q 时#!!7.0^)".#!7#0^I ""7收敛到零矩阵#也就是对任意的7!*#!!7.0^I ".#!7#0^I ""7有界1同理可证!!7.0^C ".#!7#0^C ""7和!!7.0^J".#!7#0^J""7对任意的7!*有界.这就证明了差分格式!+"是无条件稳定的.定理.%设@!I ;#C G #J 9#"7"是问题!#"^!,"的解#@7;#G #9是差分格式!+"的解1记<7;#G #9-@!I ;#C G #J 9#"7".@7;#G #9#<7-!!(<7;#G #9);-##)#-#*#.#"G -##)#-#*).#"9-##)#-#*,.##那么存在一个正常数$#使得%<7%K$!.)#3)##3))#3),"1证明!记*7-!!(+7;#G #9);-##)#-#*#.#"G -##)#-#*).#"9-##)#-#*,.##那么方程!D "减去!+"的矩阵形式为!7.0^I "!7.0^C "!7.0^J "<7##-!7#0^I "!7#0^C "!7#0^J"<7#*7##1!)*"因为0^I ’0^C ’0^J 乘积可交换#根据引理#*#方程!)*"可以写为<7##-!7.0^J ".#!7#0^J "!7.0^C ".#!7#0^C ",!7.0^I ".#!7#0^I "<7#!7.0^J ".#!7.0^C ".#!7.0^I".#*7##1用矩阵的)B 范数作用上式两端#并利用引理A 和+#可得%<7%K %<*%##7.#X -*%*X ##%1因<*-$#所以存在正常数$#使得%<7%K$!.)#3)##3))#3),"1-%数值结果下面#我们通过两个数值算例验证本文所提出的数值格式的稳定性和收敛阶#也就是说格式是有效的#并在时间和空间方向都具有较高的二阶收敛精度.设"和"3分别表示问题!#"i !,"的解析解和采用格式!#*"i !#,"得到的数值解#用离散的Q q 和Q )范数计算全局截断误差#即%"3."%Q q o -I 9‘#’9’*,.##’G ’*).##’;’*#.#@,;#G #9.@!I ;#C G #J 9#Z"#%"3."%Q )o -!#*,.#9-##*).#G -##*#.#;-#@,;#G #9.@!I ;#C G #J 9#Z ")3#3)3,"#L )1!!算例&%在问题!#"i !,"中#取*-!*##"/!*##"/!*##"#Z-#$对流和扩散系数分别为X #-X )-X ,-.##:#-:)j P #-P )-2#-2)-#$初值取@*!I #C #J "-I ,!#.I ",C ,!#.C ",J ,!#.J ",1已知的解析解为@!I#C #J #""-1."I ,!#.I ",C ,!#.C ",J ,!#.J ",#由以上条件容易算出!!I #C #J #""1对常系数算例#取优化的步长比例,-*#-*)-*,进行测试.表#列出的数值结果表明#用格式!#*"i !#,"计算常系数问题!#"i !,"时#算法是无条件稳定的#而且在时间及空间方向都是二阶收敛的#这和理论分析的结果一致.算例+%在问题!#"i !,"中#取*-!*##"/!*##"/!*##"#Z -#$对流和扩散系数分别为X #-*’)"I #X )-*’)"C #X ,-*’)"J #:#-I "#:)-!#.I ""#P #-C !#P )-!#.C "!#2#-J %#2)-!#.J "%$初值取@*!I #C #J "-I ,!#.I ",C ,!#.C ",J ,!#.J ",1已知的解析解为@!I#C #J #""-1."I ,!#.I ",C ,!#.C ",J ,,!#.J ",#由以上条件!!I #C #J #""容易算出.表)列出了变系数算例)的数值结果#这里也取优化的步长比例,-*#-*)-*,进行测试#数值结果表明用格式!#*"i !#,"计算变系数问题!#"i !,"时#算法是无条件稳定的#而且在时间及空间方向也都具AE Copyright©博看网 . All Rights Reserved.!第#期聂玉峰#等$求解三维空间分数阶对流扩散方程的d%53&9S B V522格式!!表&!算例#在时刻"-##取,-*#-*)-*,的数值误差和收敛阶345*&!P011@@%@S92><%281@312<1%@>1@S\%@1‘9I L&1#9:"j#[7:0,-*#-*)-*,.%"3."%Q q+:"<%"3."%Q)+:"<"-#’) !-#’) %-#’)#Q#*#Q)*#Q E*#Q A*#Q#C*#’#*)D1X**D)’CD*#1X**AC’C,+"1X**+#’CCC,1X**+E’#D*"1X*#*X)’*EC#)’**D D#’++E E#’++A E#’+,C E1X**AE’D)E)1X**+#’#DA E1X**+)’+"))1X*#*D’,+C E1X*##X)’*,"))’**,)#’++D*#’++C+"-#’E !-#’" %-#’C #Q#*#Q)*#Q E*#Q A*#Q#C*#’)",E1X**D,’*#C+1X**AD’"*,)1X**+#’AAC E1X**+E’DE)"1X*#*X)’*"E D)’**D"#’++#+#’++#+)’*A+*1X**A"’*D+,1X**+#’)D*D1X**+,’)**C1X*#*A’*CE+1X*##X)’*E*##’+++*#’+A+)#’+AA C"-#’+ !-#’+ %-#’+#Q#*#Q)*#Q E*#Q A*#Q#C*)’,*"A1X**DE’D,+A1X**A#’*A""1X**A)’C)A"1X**+C’"#,D1X*#*X)’)A)E)’#)C")’*EC*)’*#)D,’"*"E1X**AD’,E,A1X**+#’D*,D1X**+E’#"A E1X*#*#’*,E"1X*#*X)’)""*)’#*D+)’*,E C)’**D#表+!算例)在时刻"-##取,-*#-*)-*,的数值误差和收敛阶345*+!P011@@%@S92><%281@312<1%@>1@S\%@1‘9I L&1)9:"j#[7:0,-*#-*)-*, .%"3."%Q q+:"<%"3."%Q)+:"<"-#’) !-#’) %-#’)#Q#*#Q)*#Q E*#Q A*#Q#C*#’*+)+1X**D)’C#E"1X**AC’EAE A1X**+#’C)#)1X**+E’*"A,1X*#*X)’*C,C)’*##E)’****#’++A#)’#*A C1X**AE’+C*C1X**+#’)),D1X**+,’*")*1X*#*D’C,**1X*##X)’*AD D)’*#+,)’**,E)’****"-#’E !-#’" %-#’C #Q#*#Q)*#Q E*#Q A*#Q#C*#’*A#D1X**D)’C#"E1X**AC’"#E+1X**+#’C,"D1X**+E’##*+1X*#*X)’*EA))’**")#’++,A#’++)E)’*,"A1X**AE’AE*+1X**+#’)**"1X**+,’**A D1X*#*D’"",C1X*##X)’*D)))’*##C#’++C E#’++,+"-#’+ !-#’+ %-#’+#Q#*#Q)*#Q E*#Q A*#Q#C*+’AE,"1X**A)’,C#,1X**A"’ADC E1X**+#’EA*D1X**+,’DEC#1X*#*X)’*"+C)’**C C#’+AA D#’+A)A#’C+*C1X**AE’*C#"1X**+#’*#,)1X**+)’""A,1X*#*C’EA,"1X*##X)’*"D")’**,##’+A"D#’+A*,有二阶收敛率#这和理论分析的结果是非常吻合的..%结论本文将求解三维整数阶抛物方程的经典d%53&9S B V522格式推广到分数阶#提出了一种求解三维;F G d=的有效的数值方法#并证明该格式具有无条件稳定性和较高的二阶收敛精度#必要而充足的数值实验验证了理论结果.最后#由于分数阶导数是非局部算子#对于多维空间问题的求解需要耗费较大的计算工作量和空间存储量#在今后的工作中#我们将考虑开展适当的快速算法#以减少计算花费和加快计算速度.参考文献!(#)!;e6;G=.G21‘L&7<7:0730%@>1@I1:0%>\%@\@9<:7%29&9>81<:7%2>7\\5S7%21O59:7%2S(-).-%5@29&%\<%I L5:9:7%29&L0M S7<S# )*#E#)DA&)"D X)DE.())!/T G(VT#]^6F#/T6G(V U#1:9&.(5I1@7<9&929&M S7S%\921[S L9<1B:7I189@79Z&1\@9<:7%29&%@>1@9>81<:7%2B>7S L1@S7%2 1O59:7%2(-).G LL&71>I9:01I9:7<S92><%I L5:9:7%2#)*#E#)E)&"E#X""*.+E Copyright©博看网 . All Rights Reserved.*"郑州大学学报!理学版"第"#卷(,)!=_$^($-#T=6=_(#_e e U-U.(5I1@7<9&9LL@%‘7I9:7%2%\9:7I1>1L12>12:#2%2&7219@#S L9<1B\@9<:7%29&>7\\5S7%21O59:7%2 (-).;^G?R%5@29&%225I1@7<9&929&M S7S#)**D#E"!)"&"D)X"+#.(E)!/T=(Vff#]^YU#/T G e/V.G2%:1%2:01\727:11&1I12:I1:0%>\%@:01S L9<1\@9<:7%29&9>81<:7%2>7\\5S7%21O59:7%2(-).Y%I L5:1@S92>I9:01I9:7<S[7:09LL&7<9:7%2S#)*#*#"+&#D#A X#D)C.(")!T=-G/^T#?e_e(=fP#]^6F.;:9Z7&7:M92><%281@312<1%\9\727:18%&5I1I1:0%>\%@:01S L9<1\@9<:7%29&9>81<:7%2B>7S L1@B S7%21O59:7%2(-).-%5@29&%\<%I L5:9:7%29&92>9LL&71>I9:01I9:7<S#)*#E#)""&CAE X C+D.(C)!虎晓燕#韩惠丽.重心插值配点法求解分数阶F@1>0%&I积分方程(-).郑州大学学报!理学版"#)*#D#E+!#"&#D X),.(D)!/T=(V?#]^6F#G(T$#1:9&.G0730B%@>1@S L1<:@9&I1:0%>\%@:01I5&:7B:1@I:7I1B\@9<:7%29&>7\\5S7%21O59:7%2S(-).G LBL&71>I9:01I9:7<9&I%>1&&723#)*#C#E*!D Q A"&E+D*X E+A".(A)!d=(Va T#Y T=(?T.=\\7<712:25I1@7<9&9&3%@7:0I S\%@:0@11B>7I12S7%29&\@9<:7%29&L9@:79&>7\\1@12:79&1O59:7%2S(-).-%5@29&%\<%I L5:9:7%29&I9:01I9:7<S#)*#E#,)!E"&,D#X,+#.(+)!Y T=(-#]^6F#]^6h#1:9&.(5I1@7<9&S7I5&9:7%2\%@:01:0@11B>7I12S7%2\@9<:7%29&S5ZB>7\\5S7%21O59:7%2(-).G LL&71> I9:01I9:7<9&I%>1&&723#)*#E#,A!#""&,C+"X,D*".(#*)aG(VT#d6(.F9S:9&:1@29:723B>7@1<:7%2\727:1>7\\1@12<1I1:0%>S\%@:0@11B>7I12S7%29&S L9<1B\@9<:7%29&>7\\5S7%21O59:7%2S (-).-%5@29&%\<%I L5:9:7%29&L0M S7<S#)*#E#)"A&,*"X,#A.(##)Y T=(?T#d=(Va T.GS1<%2>B%@>1@25I1@7<9&I1:0%>\%@:[%B>7I12S7%29&:[%B S7>1>S L9<1\@9<:7%29&<%281<:7%2>7\\5S7%2 1O59:7%2(-).G LL&71>I9:01I9:7<9&I%>1&&723#)*#E#,A!#,"&,)EE X,)"+.(#))]G6HG-.?9:@7‘929&M S7S\%@S<712:7S:S92>1237211@S(?).U07&9>1&L079&;^G?#)**".(#,)h6G_P=_e(^G#;G Y Y e_#;G]=_^F.(5I1@7<9&I9:01I9:7<S(?).H1@&72&;L@7231@B$1@&93#)**D.!9D)?4>=U D::’(:(;<!(H H<B<:E<F E N<@<H9B3N B<<=C(@<:>(9:4?F T4E<’B4E;(9:4?7C Q<E;(9:!(H H D>(9:L X D4;(9:(^=f5\123##T6-79057##)#aG(V-523923#!#’+<B<:823$<7"<8!48$49T@":";47:5)2;<72<#*48"3E<B"<87(45C"<237;2:5&7;A<8B;"C#M;N:7D#*#)+#$3;7:$ )’$455<?<4!)2;<72<#\<7:7&7;A<8B;"C4!Z<237454?C#O3<7?J34@E"***##$3;7:"75>;B4E;&d51:%:012%2B&%<9&7:M%\\@9<:7%29&>1@789:781S#\@9<:7%29&L9@:79&>7\\1@12:79&1O59:7%2S[1@1Z1::1@:%>1S<@7Z192%I9&%5S>7\\5S7%2L012%I129:092%:01@I1:0%>S.T%[181@#[07&112R%M723:01<%2B 812712<1\@%II9:01I9:7<9&I%>1&723#7:9&S%<95S1>&%:S%\:@%5Z&11S L1<79&&M72S%&8723I5&:7>7I12S7%29& <9S1S.G21\\7<712:25I1@7<9&9&3%@7:0I[9S L@%L%S1>\%@S%&8723:01:0@11B>7I12S7%29&S L9<1\@9<:7%29&9>81<:7%2>7\\5S7%21O59:7%2!;F G d="ZM3121@9&74723:01d%53&9S B V522S<01I1.;:9Z7&7:M92><%281@B312<1%\:01I1:0%>[1@1L@%81>ZM:01I9:@7‘I1:0%>.P01>1@781>9&:1@29:723>7@1<:7%27I L&7<7:!G d^"\727:1>7\\1@12<1S<01I109>:01S1<%2>%@>1@9<<5@9<M72Z%:0:7I192>S L9<1>7@1<:7%2S#@1S L1<:781&M.P011\\7<712<M92><%281@312<1%@>1@S[1@1\729&&M>1I%2S:@9:1>ZM S%I125I1@7<9&1‘9I L&1S.J<G K9B C>&:0@11B>7I12S7%29&;F G d=$G d^S<01I1$Y@92JB(7<%&S%2S<01I1$d%53&9S B V522S<01I1$ S:9Z7&7:M$<%281@312<1!责任编辑&方惠敏"Copyright©博看网 . All Rights Reserved.。
基于分数阶导数的岩石非线性蠕变损伤模型

基于分数阶导数的岩石非线性蠕变损伤模型王晴;仇晶晶;朱其志;刘思利;余健【摘要】利用Riemann-Liouville分数阶理论,给出一种分数阶软体元件及其本构方程,阶数取值不同,可分别模拟蠕变的三个阶段.采用两个分数阶软体元件与虎克体进行组合,引入岩石硬化函数、损伤变量,提出一种新的含分数阶导数的非线性蠕变损伤模型,并推导出该模型的本构方程.利用砂岩的蠕变试验数据进行验证,发现该模型能有效描述砂岩的蠕变特性.【期刊名称】《河南科学》【年(卷),期】2019(037)003【总页数】5页(P406-410)【关键词】分数阶;软体元件;非线性蠕变【作者】王晴;仇晶晶;朱其志;刘思利;余健【作者单位】河海大学岩土力学与堤坝工程工程教育部重点实验室,南京 210098;河海大学江苏省岩土工程技术工程研究中心,南京 210098;江苏纬信工程咨询有限公司,南京 210014;河海大学岩土力学与堤坝工程工程教育部重点实验室,南京210098;河海大学江苏省岩土工程技术工程研究中心,南京 210098;河海大学岩土力学与堤坝工程工程教育部重点实验室,南京 210098;河海大学江苏省岩土工程技术工程研究中心,南京 210098;河海大学岩土力学与堤坝工程工程教育部重点实验室,南京 210098;河海大学江苏省岩土工程技术工程研究中心,南京 210098【正文语种】中文【中图分类】TU45在工程实际中,岩石蠕变特性是岩石材料重要力学性质之一,与岩体工程的安全性、稳定性密切相关,因此,岩石蠕变研究越来越被人们重视[1].国内外研究者多采用元件模型来研究蠕变本构关系,代表性的有Burgers模型、Maxwell模型等,其中的胡克元件仅可描述纯弹性行为,可模拟蠕变线弹性阶段;牛顿体元件仅可描述材料的黏滞性,用以模拟蠕变稳定阶段.基于此,国内外很多学者引进非线性元件来模拟蠕变加速阶段,并适当组合非线性元件和线性元件,从而能更好地模拟蠕变的三个阶段[2-11].为更好地模拟,必须使用多个元件,导致了模型参数增加.由此,一些学者提出了含参数少、方程简洁的分数阶软体元件.殷德顺等[12]基于分数阶微积分理论,提出一种介于理想固体和流体之间的分数阶软体元件,将该软体元件与弹簧元件分别进行串联和并联得到两种模型,推导出本构方程,与土的试验数据进行比较,发现模型能有效描述蠕变特性;何志磊等[13]改进了西原模型,将西原模型中的牛顿黏壶替换为分数阶Abel黏壶,克服了西元模型不能描述蠕变加速阶段的缺点,得到分数阶非定常蠕变模型;吴斐等[14]引入了分数阶非线性黏壶元件,提出了参数少且可描述蠕变三个阶段的新的岩石蠕变模型;郭佳奇[15]等将Kelvin-Voigt蠕变模型中的牛顿体元件替换为分数阶软体元件,得到分数阶Kelvin-Voigt蠕变模型,与实验数据拟合后,比较发现分数阶Kelvin-Voigt蠕变模型不仅参数更少,且拟合效果比整数阶5参数开尔文蠕变模型和整数解Kelvin-Voigt蠕变模型好.在上述研究基础上,本文讨论了阶数取值对分数阶软体元件的影响,提出了一种新的含有分数阶导数的非线性蠕变损伤模型,利用两个不同阶数的分数阶软体元件,同时引入岩石硬化函数和损伤变量,最后结合砂岩蠕变试验进行验证.1 分数阶非线性蠕变损伤模型建立1.1 含分数阶导数的非线性蠕变损伤模型分数阶微积分定义有多种,本文采用Riemann-Liouville型分数阶微积分算子理论,Riemann-Liouville积分形式为式中:为Gamma函数[15-16].Riemann-Liouville微分形式为函数 f(t)的r阶Riemann-Liouville分数阶微分[15-16].将应力-应变关系用分数阶微分形式进行表示,如式(3):式中:ξ和β均为材料常数.可知,当β=0,式(3)可表示理想固体的应力-应变关系:σ(t)-ε(t),满足胡克定律;当β=1,式(3)可表示理想流体的应力-应变关系:σ(t)-d1ε(t)/dt1,满足牛顿黏性定律.当σ(t)=const应力恒定时,可描述蠕变现象.对式(3)两边进行积分,可得软体元件的蠕变方程为[17-18]:如图1和图2所示,应力水平保持恒定,当β采用不同值时,由式(4)表示的一系列蠕变曲线.可知,随着β取值的增大,软体元件应力-应变的非线性增加,这与理想固体保持不变的状态和理想流体线性增加的状态不同.从图1可以看出,当0<β<1时,随着β增加,软体元件的应力应变关系线性特征逾明显,可以用来描述蠕变的线性、稳态阶段;从图2可以看出,当β>1时,随着β增加,软体元件的应力应变曲线的斜率越来越大,可以描述蠕变加速阶段.图1 分数阶软体元件的蠕变曲线图(0<β<1)Fig.1 Creep curves of fractionalsoftware components(0<β<1)由此发现,可以通过调整分数阶的阶数来模拟蠕变的不同阶段,故提出如图3所示模型:第I部分为引入硬化函数的虎克体元件;第II部分为黏弹性分数阶软体元件,其阶数小于1;第III部分为引入损伤变量、阶数大于1的黏塑性分数阶软体元件和应力阈值开关元件并联组成.图2 分数阶软体元件的蠕变曲线图(β>1)Fig.2 Creep curves of fractional software components(β>1)图3 含分数阶导数的非线性蠕变损伤模型Fig.3 Nonlinear creep damage model with fractional derivatives1.2 蠕变本构方程对于第I部分,引入硬化函数,由陈化理论[19]可知,材料强度H(σ,t)可表示为式中:H(σ ,t)为材料强度;σ、t分别为应力、时间;H0为初始强度;λ、n为材料常数.本文为方便计算,取n=0.则第I部分的本构方程为对于第II部分,使用黏弹性分数阶软体元件,其阶数小于1,本构方程为式中:γ1,α分别为材料常数和分数阶阶数,且0<α<1.对式(7)两边积分,可得第II部分本构方程为对于第III部分,引入损伤变量.根据能量损伤的方法,定义损伤变量为:式中:D(σ,t)表示t时刻材料的损伤变量;E(σ,t)表示任意t时刻的弹性模量,其值与t时刻对应的应力水平有关;E0表示材料初始弹性模量.由经典蠕变损伤模型可定义E(σ,t)为式中:e为自然常数;θ表示材料常数.将式(10)代入式(9),可得由式(11)可知:t→∞时,D=1,表示材料已经完全损伤.易知,随时间t和应力σ增长,试样弹性模量逐渐减小,损伤慢慢变大.Kachanov定义有效应力为式中:为有效应力;σ为名义应力.将式(11)带入式(12),得当σ>σs(长期强度)时,进入蠕变损伤阶段,则本构方程为式中:γ2,β分别为软体元件的材料常数和阶数,且β>1.对式(14)积分,并将式(13)代入,得第III部分本构方程为综上,当σ<σs时,蠕变无加速阶段,模型中仅存在第I、II部分,此时本构方程为当σ>σs时,蠕变会进入加速阶段,此时模型包含第I、II、III部分,本构方程为2 蠕变实验及模型验证2.1 实验准备及实验结果试验采用全自动岩石三轴伺服仪,如图4所示.试样为砂岩(如图5所示),标准圆柱形,尺寸为50 mm×100 mm(直径×高度).图4 全自动岩石三轴伺服仪Fig.4 Automatic three-axis rock servo试验室温度控制在(21±0.5)℃内,以减小温度对试验数据的影响.试验期间,先将围压加载至目标值10 MPa,并保持恒定,加载速率3 MPa/min.待围压、温度稳定后,采用轴向应力控制方法加载偏压至目标值(分别为:118、119 MPa),加载速率为0.3 MPa/min.保持温度、偏压稳定,直到试样破坏.蠕变试验曲线如图6所示.图5 砂岩试样Fig.5 Sandstone sample2.2 模型验证及参数识别利用以上蠕变试验数据验证模型,对实验数据进行非线性拟合,得模型参数见表1.由于试样之间的差异性,所得参数存在一定差异.如图7所示,根据拟合结果可知,本文提出的非线性蠕变损伤模型能较好地模拟蠕变试验的全过程.图6 不同强度下的蠕变试验应力-应变曲线Fig.6 Stress-strain curve of creep test under different strength图7 非线性蠕变损伤模型的拟合曲线Fig.7 Fitting curve of nonlinear creep damage model表1 模型参数Tab.1 Parameters of model应力水平/MPa 118 119 λ α β θ0.484 0.456 γ1/(h·MPa-1)0.726×103 0.789×103 0.245 0.286 γ2/(h·MPa-1)3.967×1011 2.177×1011 8.121 9.152 0.301 0.6713 结论1)基于Riemann-Liouville分数阶微积分理论,给出了分数阶软体元件及其本构关系,并发现当分数阶阶数取值小于1时,软体元件可模拟瞬时弹性和稳定蠕变阶段,当阶数取值大于1时,软体元件可以模拟加速蠕变阶段.2)通过两个软体元件与虎克体进行组合,引入硬化函数和损伤变量,增加应力阈值开关,提出含分数阶导数的非线性蠕变损伤模型,并推导出该模型的本构方程.由砂岩试验数据进行验证,提出的非线性蠕变损伤模型能较好地描述蠕变全过程.【相关文献】[1]徐平,杨挺青.岩石流变试验与本构模型辨识[J].岩石力学与工程学报,2001,20(S1):1739-1744.[2]徐卫亚,杨圣奇,杨松林,等.绿片岩三轴流变力学特性的研究(I):试验结果[J].岩土力学,2005,26(4):531-537.[3]徐卫亚,杨圣奇,谢守益,等.绿片岩三轴流变力学特性的研究(II):模型分析[J].岩土力学,2005,26(5):23-28.[4]杨广雨,王伟,熊德发,等.岩石非线性黏弹塑性蠕变模型研究[J].河北工程大学学报(自然科学版),2017,34(4):23-26.[5]蒋昱州,张明鸣,李良权.岩石非线性黏弹塑性蠕变模型研究及其参数识别[J].岩石力学与工程学报,2008,27(4):832-839.[6]赵延林,曹平,文有道,等.岩石弹黏塑性流变试验和非线性流变模型研究[J].岩石力学与工程学报,2008,27(3):477-486.[7]张英.岩石流变的一种非线性黏弹塑性流变模型研究[J].湖南工业大学学报,2015,29(3):10-14.[8]蒋海飞,胡斌,刘强,等.一种新的岩石黏弹塑性流变模型[J].长江科学院院报,2014,31(7):44-48.[9]马明军,钟时猷.一个软弱岩石的粘弹塑性流变力学模型[J].中南矿冶学院学报,1990,21(3):236-241.[10] LIU L,WANG G,CHEN J,et al.Creep experiment and rheological model of deep saturated rock[J].Transactions of Nonferrous Metals Society of China,2013,23(2):478-483.[11] MARANINI E,YAMAGUCHI T.A non-associated viscoplastic model for the behaviour of granite in triaxial compression[J].Mechanics of Materials,2001,33(5):283-293.[12]殷德顺,任俊娟,和成亮,等.一种新的岩土流变模型元件[J].岩石力学与工程学报,2007,26(9):1899-1903.[13]何志磊,朱珍德,朱明礼,等.基于分数阶导数的非定常蠕变本构模型研究[J].岩土力学,2016,37(3):737-744.[14]吴斐,刘建锋,武志德,等.盐岩的分数阶非线性蠕变本构模型[J].岩土力学,2014,35(S2):162-167.[15]郭佳奇,乔春生,徐冲,等.基于分数阶微积分的Kelvin-Voigt流变模型[J].中国铁道科学,2009,30(4):1-6.[16]何明明,李宁,陈蕴生,等.基于分数阶微积分岩石的动态变形行为研究[J].岩土工程学报,2015,37(S1):178-184.[17]宋勇军,雷胜友.基于分数阶微积分的岩石非线性蠕变损伤力学模型[J].地下空间与工程学报,2013,9(1):91-95.[18]王学彬.拉普拉斯变换方法解分数阶微分方程[J].西南师范大学学报(自然科学版),2016,41(7):7-12.[19]穆霞英.蠕变力学[M].西安:西安交通大学出版社,1990.。
极限热弹性材料微结构的拓扑优化设计

。 明三相材料可以获得超出所组成单相材料 的特殊性能1 2 1 , 这与它 用许多更成熟 的基于梯度的优化算法 们的微结构的布局 和它们的材料组合物相关 。
到 目前为止 ,国内外 已有一些研究通过拓扑优化方法实现 拓扑优化方法是—个迭代 的数值计算过程 ,它通过在一个固 具有特定性能参数的周期性复合材料 的微结构设计 。文献【 J 较早 定的参考范围内重新分配—个给定的量的材料 , 直到有关的结构性 提出材料微结构拓扑优化理论 , 对各项同性材料微结构单胞采用
Me t h o d
1 引言
具有特殊热弹性性能的超材料一直技术研究与基础研究领
能约束最小化, 以确定最佳的拓扑结构日 。在过去的二十年中, 拓扑 优化学科经历了长足的进步, 已开发了几个不同的方法 。水平集
一
不同的材料界面形状的保真 的关注点『 l 1 。它通常通过结合三种常规材料 , 如金属或塑料 , 形成 方法 由于其光滑等设计 的边界特征 , 以及集成 的形状优化 、 拓扑优化 , 已经成 个新的周期性微结构的复合材料来实现。在工程应用中 , 温度 度和拓扑结构的灵活性,
Ma t er i a l s wi t h E x t r e me Th e r mo e l a s t i c P r o p e r t i e s
W ANG Yu ,L V E n - l i ,GUO J i a - mi n g , Z HAO J u n - h o n g
( C o l l e g e o f E n g i n e e r i n g , S o u t h C h i n a A g r i c u l t u r a l U n i v e r s i t y , G u a n g d o n g G u a n g z h o u 5 1 0 6 4 2 , C h i n a )
基于分数阶微积分理论的粘弹性流体流动与传热研究

结论
本次演示通过概述分数阶微积分理论的基本概念和其在粘弹性流体流动与传 热研究中的应用,探讨了粘弹性流体流动与传热之间的关系。利用数值模拟方法 对两者之间的关系进行计算和分析,揭示了不同参数对流体流动和传热的影响。 结果表明,分数阶微积分理论在粘弹性流体流动与传热研究中具有重要的应用价 值和广阔的发展前景。
传热对流体流动的影响:当流体的温度分布不均匀时,会导致温差应力,从 而影响流体的流动特性。分数阶微积分理论可以描述这种由于温度场引起的非线 性应力场,从而预测流体在传热作用下的流动行为。
流体流动对传热的影响:流体的流动状态会影响传热的效果。在粘弹性流体 中,由于流体的粘弹性,流体的流动性受到限制,从而影响热量的传递。分数阶 微积分理论可以描述这种由于流体流动性变化引起的传热行为的改变。
分数阶微积分理论概述
分数阶微积分理论是一种扩展经典微积分理论的方法,它允许我们将传统整 数阶微分运算的阶数扩展到任意实数,从而更好地描述各种复杂现象。分数阶微 分的定义基于函数在某一点的泰勒级数展开,根据展开式的系数来定义分数阶微 分。分数阶微积分具有一些独特的性质,如非局部性、记忆效应等,使其在处理 某些非线性问题时具有优势。
两者之间的相互作用:粘弹性流体流动与传热之间存在相互影响的关系。一 方面,传热过程会改变流体的温度分布和流动性;另一方面,流体的流动状态会 影响传热的效果。分数阶微积分理论可以描述这种复杂的相互作用过程,从而预 测和理解粘弹性流体在流动和传热过程中的整体行为。
数值模拟分析
为了深入理解和研究粘弹性流体流动与传热之间的关系,数值模拟是一种有 效的方法。通过建立数学模型,利用计算机进行数值计算和分析,可以揭示粘弹 性流体流动与传热的内在机制。
基于分数阶微积分理论的粘弹性流 体流动与传热研究
化学动力学的分数阶微积分学

化学动力学的分数阶微积分学化学动力学是研究化学反应速率与反应机理的学科。
分数阶微积分学是在传统微积分学的基础上,引入了分数阶概念,并将其应用到不同的学科领域中,包括化学动力学。
本文将探讨分数阶微积分学在化学动力学中的应用。
1. 分数阶微积分学的基本概念传统的微积分学是基于整数阶的概念,例如导数和积分。
分数阶微积分学则是引入了分数阶的概念,根据分数阶的不同,可以得到不同的导数和积分。
例如,分数阶求导可以用分数阶微分方程表示,而分数阶积分可以用分数阶积分方程表示。
分数阶微积分学的应用十分广泛,包括物理学、生物学、金融学、信号处理等领域。
在化学动力学中,分数阶微积分学的应用也得到了广泛关注。
2. 分数阶动力学方程化学反应的速率通常用速率常数表示,速率常数可以在实验中通过测量反应物消耗的速度来确定。
然而,有些反应速率并不是简单的一阶动力学反应,而是更复杂的分数阶动力学反应。
分数阶动力学反应通常由下列方程描述:$$\frac{d^{\alpha} [A]}{dt^{\alpha}} = k [A]^{\beta}$$其中,$[A]$为反应物的浓度,$k$为速率常数,$\alpha$和$\beta$为实验中测定的反应动力学指数。
当$\alpha=1$,$\beta=1$时,上述方程即为经典的一阶反应动力学方程。
然而,当$\alpha$和$\beta$分别等于$0.5$时,方程的积分形式为:$$[A](t) = A_0 \left(1 + k_1t^{0.5}\right)^{-2}$$其中,$A_0$为初始浓度,$k_1=\frac{k}{\sqrt{\pi}}\Gamma \left(\frac{3}{2}-\alpha\right)$为分数阶速率常数,$\Gamma$为伽玛函数。
3. 分数阶反应的动力学模型分数阶反应通常具有不同于一阶反应的反应动力学特征。
因此,为了更好地描述分数阶反应的动力学模型,需要引入新的数学工具。
分数阶微积分流变模型在岩体结构加速流变破坏分析中的应用

r
i《 烈 蜒蜥
鲤 霹 孵蜘彝 矗 ℃ 零 卿 辩 够 徊酗
誓
分数 阶微 积分 流 变模 型在 岩体 结构 加 速 流变 破 坏 分 析 中 的应 用
黄耀 英 , 郑 宏
40 7 ) 3 0 1
( .三 峡 大 学 水利 与 环 境 学 院 , 北 宜 昌 4 30 ; 1 湖 4 0 2
关键 词 : 数阶微积 分 ;流 变模 型 ; 速流 变 ; 分 加 破坏
中图分 类号 : U 5 T 15 T 4 7; B 1
文献标 志码 : A
Ap lc to ff a to a r e a c l s r e l g c lm o e n p ia i n o r c i n lo d r c lu u h o o i a d li
c lua in f r l a e e u e ac lto 0 mu a r d d c d. Th n r cina o d r a c l s h oo ia mo l s p le t t e e fa to l r e c lu u r e lgc l de i a p id o h
C ieeA a e yo c n e , h n4 0 7 ,C ia hn s cd m f i c s Wu a 3 0 1 h ) Se n
Absr c :Ast h ia v n a e o e d n n r mee s wh n c n e to a l me tmo li s d ta t o t e d s d a t g fn e i g ma y paa t r e o v ni n lee n de s u e t te p rme t aa wel h h o y o a t n lo d rc lu u sa p id t h h o o y a ay i f o f x e i n a d t l,t e te r ff ci a r e ac l s i p l o t e r e lg n l ss o i l r o e
利用分数阶(G′G)展式法构造分数阶KdV-Burger方程方程的精确行波解

利用分数阶(G′G)展式法构造分数阶KdV-Burger方程方程的精确行波解尹伟石;李琰;徐飞【摘要】(G′G)展式法是一种行之有效的求解分数阶偏微分方程的方法.利用行波变化与齐次平衡技巧可以对该方法进行拓展,拓展后的方法能够处理更一般的分数阶偏微分方程.最后将拓展后的方法应用到基于黎曼-刘维尔积分意义下的时间空间分数阶KdV-Burger方程中,通过符号计算可以得到方程的精确行波解.与其他方法相比,拓展的(G′G)展式法不需要进行变换和数值逼近,计算更加的简洁.%(G′ G) expansionmethod is an effective method for solving fractional partial differentialequations.The method can be extended by using the traveling wave variation andthe homogeneous balance technique,and the extended method can be used to dealwith the more general fractional partial differential equations.Finally, theextended method is applied to the time space fractional KdV-Burger equationbased on the Liu Weier Riemann integral, and the exact traveling wave solutionsof the equations can be obtained by the symbo lic computation. Compared withother methods, (G′ G) ex-pansionmethod don't need to doing transform and numerical approximation,so thecalculation is more simple.【期刊名称】《长春理工大学学报(自然科学版)》【年(卷),期】2016(039)006【总页数】4页(P125-128)【关键词】分数阶(G′G)展式法;分数阶KdV-Burger方程;精确行波解【作者】尹伟石;李琰;徐飞【作者单位】长春理工大学理学院,长春 130022;长春理工大学理学院,长春130022;东北师范大学数学与统计学院,长春 130024【正文语种】中文【中图分类】O241.82近年来,分数阶偏微分方程(FPDEs)频繁地出现于物理、生物、工程、信号处理、系统识别、控制理论、金融和分子动力学等领域,已经成为偏微分方程领域关注的焦点问题。
基于分数阶微积分的KelvinVoigt流变模型

1、KEELVIN-VOIGT流变模型
KelvinVoigt流变模型是一种描述材料在应力作用下的变形行为的本构方程, 由Kelvin和Voigt于20世纪初提出。该模型基于弹簧和粘壶串联的物理模型,将 材料的弹性行为和粘性行为相结合,适用于描述多种材料的流变特性,如橡胶、 混凝土、生物组织等。
在KelvinVoigt流变模型中,材料在受到外部应力作用时,其变形行为可以 描述为:
此外,分数阶微积分理论在KelvinVoigt流变模型中的应用还具有以下优越 性:
(1)能够更加精确地描述材料的蠕变和松弛行为; (2)可以根据实验数据确 定模型参数,提高模型的预测精度; (3)可以建立材料在不同应力水平下的响应 模型,更好地理解材料的力学行为; (4)可以在一定程度上克服经典Kelvin模型 和Voigt模型的局限性,扩展其应用范围。
3、CASE STUDY:KelvinVoigt 流变模型的应用
为了验证KelvinVoigt流变模型的有效性和适用性,本节通过案例分析对其 进行探讨。
首先,在材料科学领域,KelvinVoigt流变模型被广泛应用于橡胶、塑料等 高分子材料的研究。例如,在对橡胶进行疲劳试验时,该模型可以准确描述其在 不同应力水平下的疲劳行为,为橡胶制品的优化设计和疲劳寿命预测提供了有力 支持。
变参数黄土蠕变损伤模型
基于分数阶微积分的变参数黄土蠕变损伤模型,主要考虑了黄土的蠕变特性、 损伤特性和环境因素。这个模型采用分数阶微积分描述黄土的蠕变行为,同时引 入了损伤变量来描述黄土的损伤演化。在模型中,我们采用了可变的参数来描述 黄土在不同环境条件下的蠕变和损伤特性。这些参数可以通过实验和数据分析得 到,使得模型更加灵活和实用。
该模型的优点在于: 1、更好地描述了黄土的蠕变和损伤特性; 2、考虑了环境因素的影响;
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
为求解式 ( 18) 具有复杂支点像函数的反演, 应用文献 [ 11] 提出的离散化求 L aplace 逆变换的
( 8)
式中 : f ( s) 为 f ( t ) 的 L aplace 变换。
当R ( t) = R 0H ( t) ( R 0 为初始恒定应力, H ( t) 为 H eaviside 单位阶跃函数 ) 时, 式 ( 8) 经 Laplace 逆变换后, 整理得到 F C 元件的蠕变方程:
0 tB E ( t) = R N #( 1 + B )
( t) 。
当 f ( t ) 在 t = 0 附近可积, 0 [ B [ 1, 分数阶 微积分的 L aplace 变换公式为
B -B L [ 0 Dt f ( t ) , s] = s f ( s) B L [ 0 DB t f ( t) , s] = s f ( s)
B> 0 B> 0
( 4)
岩土的流变特性是岩土介质的重要力学性质之 一, 岩土工程的长期稳定性和耐久性都与岩土介质 的流变性有密切关系 [ 1, 2] 。例如, 深埋地下的洞室 或巷道在成洞之初, 周围岩体呈现稳定 , 随时间推 移, 岩体的变形不断发展 , 经过一些时日, 洞体可 能失稳或坍塌破坏, 围岩具有随时间增长而缓慢变 形的明显特征; 软土地基结构物内力随时间增长变 化的长期、复杂过程 , 也是软土时效流变特性的表 现 。在工程实践中, 岩土的流变现象包括蠕变、 应力松弛、流动、长期强度等。随着人们对岩土工 程长期安全的重视 , 岩土的流 变研究越来越 被重 视。 岩土流变本构模型的研究就是探讨用何种本构 方程描述岩土材料的应力、应变和时间的关系。方 程既要能较准确地反映岩土的流变特性 , 又要符合 工程实际。目前使用的流变模型主要有经验流变模 型、组合流变模型、内时模型和蠕变损伤断裂模型 等。其中组合流变模型是研究流变问题最常使用的 模型 , 它按照岩土的弹性、塑性和黏滞性设定一些 基本元件 , 然后根据岩土的具体性质, 通过调整模 型的参数和组合元件的数目, 将其组合成能基本反
1
分数阶微积分理论
分数阶微积分是 1 个研究任意阶次的微分、积
分算子 特 性及 应 用的 数 学 问题。本 文 采 用 Rie mann - L io uvill e 型分数阶微积分算子理论 , 定义函 数的任意阶微分和积分公式 [ 7] 。
2
中
国
铁
道
科
学
第 30 卷
对于任意复数 B, Re( B ) > 0, 函数 f ( t ) 的 B阶 积分定义为 d f ( t) = -B dt 其中 , #( B ) =
3
3 11
分数阶流变模型
经典整数阶 Kelvin -Voigt 流变模型 经典整数阶 Kelvin - Vo ig t 模型是由 1 个虎克体
元件 ( 弹性模量为 E 2 ) 与 1 个牛顿体元件 ( 黏滞 系数为 G) 并联后再与另一虎克体元件 ( 弹性模量 为 E 1 ) 串联组成 , 如图 4 所示。
图5
分数阶 K elvin - V o igt 模型
3 12 11
分数阶 Kelvin - Voigt 模型的本构方程
各元件串联时 , 各元件应力相等 , 应变等于各 元件应变之和 , 各元件并联时, 各元件应变相等, 应力等于各元件应变之和[ 10] , 基于上述原则, 有 如下关系式 : 由上述分 析知, F C 元件既能 控制变形速率 , 也能控制应力 ( 应变 ) , 而牛顿体元件只能通过黏 滞系数控制变形速率, 从这一点上看 , FC 元件将 能更好地反映流变问题的非线性渐变过程。 E= E 1 + E 2 R = E1 E 1 R = E2 E 2 + N 0D E
[ 5] B B
式中 : R ( t) , E ( t) 和 E 分别为 理想弹性 体的应力、
式中 : G 为黏滞系数。 如果将式 ( 5) 改写为 R( t) = E d0 E ( t) / dt0 则黏弹性体的应力 ) 应变关系可以表述为 R( t) = N dE ( t) / dt
B B
整 FC 元件的参数 B, N改变其蠕变曲线或松弛曲线 的线型 , 从而精确拟合材料的实验结果, 使它更真 ( 7) 实地体现材料的性质。与 F C 元件比较, 牛顿体元 件的蠕变柔量与 t 成正比( J ( t) = t / G) , 松弛模量为 脉冲函数 , 从无穷大突然松弛到零( G( t) = G D ( t) ) , 牛顿体元件无法通过调整参数改变其蠕变曲线或松
1 1
蠕变柔量为 1 tB J( t) = N # ( 1 + B) 初始应变) , 可推得 F C 元件的松弛方程: R ( t) = N E 0 t # ( 1 - B)
-B
( 10)
同理, 当 E ( t) = E 0H ( t) ( E 0 对应于初始应力的
( 11)
松弛模量为 ( 5) t- B ( 12) #( 1 - B ) 不同材料参数时 F C 元件的蠕变柔量与松弛模 G( t ) = N 量随时间的变化曲线如图 2 和图 3 所示。 ( 6) 由图 2 可见, F C 元件的蠕变柔量按正分数幂 律 t 增长。由图 3 可见 , 松弛模量是从无穷大按负 分数幂律 t 松弛到零的。对于不同材料, 可通过调
B t 2
( 16)
对上述各等式分别进行 L aplace 变换, 然后整理再 进行 Laplace 逆变换 , 可得分数阶 Kelvin - Voigt 模 型的流变本构方程 : R+ 3 12 12 N B E1 E2 E1 N B 0Dt R = E+ 0Dt E E 1 + E2 E1 + E2 E1 + E2 ( 17) 分数阶模型的蠕变方程及蠕变柔量 令R ( t ) = R0 H ( t) , 代入 式 ( 17) , 两边 进行 L aplace 变换整理得 E ( t) = R 0 1 R 0 1 + E1 s E2 + N sB s ( 18)
第 3 0卷 , 第 4期 20 09 年 7月
文章编号 : 1001 - 4632 ( 2009) 04 - 0001 - 06
中 国 铁 道 科 学 CH INA RAIL WAY SCIEN CE
Vo l 1 30 No1 4
Jul y, 2009
基于分数阶微积分的 Kelvin -Voigt 流变模型
(B ) B
B
F C 元件的本 构方程定义为 式 ( 7) , 当 B= 0 时 , F C 元件就是虎克体元件; 当 B= 1 时 , FC 元 件就是牛顿体元件。 利用 Riemann - Liouville 型分 数阶微积 分算子 理论, 将公式 ( 7) 两 边同 时进 行 L aplace 变 换, 然后整理得 s- B E= R N
收稿日期 : 2008 -08 -20; 修订日期 : 2009 - 03 -20 基金项目 : 国家自然科学基金资助项目 ( 50478061) 作者简介 : 郭佳奇 ( 1981 ) ) , 男 , 河南周口拟实 际岩土的应 力 ) 应变关系。这样建立起来的模型属于一维流变 模型, 如 Max w ell 模型、 Bingham 模型、 Kelvin Voig t 模型和 Bur ger 模型等。 组合流变模型的流变本构方程是 1 种整数阶微 分形式的本构关系 , 为了很好地拟合试验结果 , 需 要较多的元件, 自然增加了本构方程中的参数。而 作为分形几何和分数维动力学基础的分数阶微积分 在建模时需要较少的参数且方程简明 [ 4] , 目前已在 松弛、随机扩散、黏弹 性力学等诸多 领域得到应 用。鉴于分数阶微积分在黏弹性本构建模方面具有 的一些优越性[ 5, 6] , 本文用含分数阶导数的力学单 元取代 Kelv in - V oigt 模型中的牛顿体元件 , 建立基 于分数阶微积分的 Kelv in - V oigt 模型。
t
0
论上能刻画参数 B无穷集合中的任何性态。 根据上述观点 , 仿照基本流变元件定义 1 种新 的力学元件 , 命名为 F C 元件 , 如图 1 所示。它能 刻画材料从理想弹性体到理想黏性体的所有性态, 含有分数阶微积分的思想。
-B
t
0
B Dt f ( t) =
Q
t
0
t
( t - S) #( B )
郭佳奇 , 乔春生 , 徐
( 1. 北 京交通大学 隧道与岩土工程研究所 , 北京 摘
1 1
冲 , 黄山秀
1
2
100044; 2. 河南理工大学 材料学院 , 河南 焦作
454000)
要 : 为研究岩土材料的应力、应变和时间 的关系 , 基 于分数 阶微积 分理论 , 定义含分 数阶导 数的力 学
元件 ( FC 元件 ) , 推导 FC 元件的蠕变柔量和松弛模量。与 牛顿体元 件相比 , F C 元件 能更好地 反映流 变问题 的 非线 性渐变过程。借鉴经典元件组合模型的建模思路 , 用 F C 元 件取代整数阶微积分 Kelvin - V oig t 流变模型中 的 牛顿 体元件 , 形成基于分数阶微积分的 K elv in - Vo ig t 流变模型。应用离散化求 Laplace 逆变换的方 法以及 H - F ox 函数 , 得出分数阶微积分 K elv in - Vo igt 流变模 型的本构 方程、蠕 变方程、松 弛方程、 蠕变柔 量及 松弛模 量的 解 析表达式。采用整数阶 微积 分 K elvin - V oig t 流变 模型、 整数 阶 5 参 数开 尔 文流 变模 型和 分数 阶微 积分 K elv in Vo ig t 流变模型对试验数据拟合的结果表明 , 分数阶微积分 Kelvin - V oigt 流 变模型不 但拟合精度 高 , 能够克服 整 数阶微积分 K elvin - V o igt 流变模型在蠕变初期及 蠕变曲 线拐点 附近与 试验数 据不能 很好吻合 的弊端 , 而且能 够 在保证拟合精度的条件下 , 减少本构模型中的参数。 关键词 : 岩土 ; 流变性 ; 蠕变 ; 分数阶微积分 ; 力学元件 ; 流变模型 中图分类号 : T U 470 文 献标识码 : A