text:
flavor 1
srclang 3
id 65535
numfuncs 5
import "RCCycleTest.mplt"
import "/home/bravewtz/Desktop/openarkcompiler/libjava-core/libjava-core.mplt"
entryfunc &LRCCycleTest_3B_7Cmain_7C_28ALjava_2Flang_2FString_3B_29V
fileinfo {
  @INFO_filename "RCCycleTest.jar"}
srcfileinfo {
  1 "RCCycleTest.java"}
type $__class_meta__ <struct {
  @shadow <* void> final,
  @monitor i32 final,
  @classloader u16 final,
  @objsize u16 final,
  @itab <* void> final,
  @vtab <* void> final,
  @gctib <* void> final,
  @classinforo <* void> final,
  @clinitbridge <* void> final}>
type $__class_meta_ro__ <struct {
  @classname <* void> final,
  @ifields <* void> final,
  @methods <* void> final,
  @superclass_or_componentclass <* void> final,
  @numoffields u16 final,
  @numofmethods u16 final,
  @flag u16 final,
  @numofsuperclasses u16 final,
  @padding u32 final,
  @mod i32 final,
  @annotation i32 final,
  @clinitAddr i32 final}>
type $__method_info__ <struct {
  @method_in_vtab_index i64 final,
  @declaringclass i64 final,
  @addr i64 final,
  @mod i32 final,
  @methodname i32 final,
  @signaturename i32 final,
  @annotationvalue i32 final,
  @flag u16 final,
  @argsize u16 final,
  @padding u32 final}>
type $__method_info_compact__ <struct {
  @method_in_vtab_index i32 final,
  @addr i32 final,
  @lebPadding0 u8 final}>
type $__field_info__ <struct {
  @offset u64 final,
  @mod u32 final,
  @flag u16 final,
  @index u16 final,
  @typeName i64 final,
  @fieldname u32 final,
  @annotation u32 final,
  @declaringclass <* void> final}>
type $__field_info_compact__ <struct {
  @offset u32 final,
  @lebPadding0 u8 final}>
type $__superclass_meta__ <struct {
  @superclassinfo <* void> final}>
type $MUIDFuncDefTabEntry <struct {
  @funcUnifiedAddr u64}>
type $MUIDFuncInfTabEntry <struct {
  @funcSize u32,
  @funcName u32}>
type $MUIDFuncDefMuidTabEntry <struct {
  @muidLow u64,
  @muidHigh u64}>
type $MUIDDataDefTabEntry <struct {
  @dataUnifiedAddr u64}>
type $MUIDDataDefMuidTabEntry <struct {
  @muidLow u64,
  @muidHigh u64}>
type $MUIDUnifiedUndefTabEntry <struct {
  @globalAddress u64}>
type $MUIDUnifiedUndefMuidTabEntry <struct {
  @muidLow u64,
  @muidHigh u64}>
type $MUIDRangeTabEntry <struct {
  @tabBegin <* void>,
  @tabEnd <* void>}>
javaclass $LRCCycleTest_3B <$LRCCycleTest_3B> public
javaclass $LTest__A_3B <$LTest__A_3B>
javaclass $LTest__B_3B <$LTest__B_3B>
var $_C_STR_d4b6a4788952b3ae6554bbc29994ed26 fstatic <[5] u64> readonly = [0, 0x2900000000, 0x73616c63ae9a5a99, 0x6320736920422073, 0x64657463656c6c6f]
var $_C_STR_a360f9874cf7de12fce9dffe2c4c052a fstatic <[6] u64> readonly = [0, 0x3100000000, 0x73616c634a44d58c, 0x6e20736920422073, 0x656c6c6f6320746f, 0x64657463]
var $Ljava_2Flang_2FSystem_3B_7Cout extern <* <$Ljava_2Fio_2FPrintStream_3B>> final public static
var $__cinf_Ljava_2Flang_2FString_3B <$__class_meta__>
func &MCC_GetOrInsertLiteral () <* <$Ljava_2Flang_2FString_3B>>
var $__vtb_LRCCycleTest_3B fstatic <[11] <* void>> = [16, 24, 32, 36, 8, 4, 48, 12, 20, 44, 28]
var $__vtb_LTest__A_3B fstatic <[12] <* void>> = [16, 24, 32, 36, 8, 4, 48, 12, 20, 44, 28, 18]
var $__vtb_LTest__B_3B fstatic <[11] <* void>> = [16, 24, 32, 36, 8, 4, 48, 12, 20, 44, 28]
var $__cinf_LRCCycleTest_3B <$__class_meta__> public
var $__methods_info__LRCCycleTest_3B fstatic <[2] <$__method_info__>> public = [[1= 0xfff6, 2= addrof ptr $__cinf_LRCCycleTest_3B, 3= addroffunc ptr &LRCCycleTest_3B_7Cmain_7C_28ALjava_2Flang_2FString_3B_29V, 4= 9, 5= 60, 6= 80, 7= 172, 8= 0xa201, 9= 1, 10= 0], [1= 0xfff6, 2= addrof ptr $__cinf_LRCCycleTest_3B, 3= addroffunc ptr &LRCCycleTest_3B_7C_3Cinit_3E_7C_28_29V, 4= 0x10001, 5= 188, 6= 216, 7= 172, 8= 0xad81, 9= 1, 10= 0]]
var $__cinf_Ljava_2Flang_2FObject_3B extern <$__class_meta__> public
var $__superclasses__LRCCycleTest_3B fstatic <[1] <$__superclass_meta__>> public = [[1= 0x4000000000000000]]
var $__classinforo__LRCCycleTest_3B fstatic <$__class_meta_ro__> public = [1= 4, 2= 0, 3= addrof ptr $__methods_info__LRCCycleTest_3B, 4= addrof ptr $__superclasses__LRCCycleTest_3B, 5= 0, 6= 2, 7= 0, 8= 1, 9= 0, 10= 33, 11= 172, 12= 0]
var $MCC_GCTIB__LRCCycleTest_3B fstatic <* void> public
var $classStateInitialized u64
var $__cinf_LRCCycleTest_3B <$__class_meta__> public = [1= 0x6f15ad0d7108f918, 2= 0, 3= 0xffff, 4= 0, 5= 0, 6= addrof ptr $__vtb_LRCCycleTest_3B, 7= addrof ptr $MCC_GCTIB__LRCCycleTest_3B, 8= addrof ptr $__classinforo__LRCCycleTest_3B, 9= addrof ptr $classStateInitialized]
var $__cinf_LTest__A_3B <$__class_meta__> public
var $__fields_info__LTest__A_3B fstatic <[1] <$__field_info__>> public = [[1= 4, 2= 0, 3= 0x6c0, 4= 0, 5= 280, 6= 268, 7= 172, 8= addrof ptr $__cinf_LTest__A_3B]]
var $__methods_info__LTest__A_3B fstatic <[2] <$__method_info__>> public = [[1= 11, 2= addrof ptr $__cinf_LTest__A_3B, 3= addroffunc ptr &LTest__A_3B_7CgenCycle_7C_28_29V, 4= 0, 5= 316, 6= 216, 7= 172, 8= 0x3f80, 9= 1, 10= 0], [1= 0xfff6, 2= addrof ptr $__cinf_LTest__A_3B, 3= addroffunc ptr &LTest__A_3B_7C_3Cinit_3E_7C_28_29V, 4= 0x10000, 5= 188, 6= 216, 7= 172, 8= 0xad81, 9= 1, 10= 0]]
var $__superclasses__LTest__A_3B fstatic <[1] <$__superclass_meta__>> public = [[1= 0x4000000000000000]]
var $__classinforo__LTest__A_3B fstatic <$__class_meta_ro__> public = [1= 232, 2= addrof ptr $__fields_info__LTest__A_3B, 3= addrof ptr $__methods_info__LTest__A_3B, 4= addrof ptr $__superclasses__LTest__A_3B, 5= 1, 6= 2, 7= 0, 8= 1, 9= 0, 10= 32, 11= 172, 12= 0]
var $MCC_GCTIB__LTest__A_3B fstatic <* void> public
var $__cinf_LTest__A_3B <$__class_meta__> public = [1= 0x13bdf92a2d6fce4b, 2= 0, 3= 0xffff, 4= 0, 5= 0, 6= addrof ptr $__vtb_LTest__A_3B, 7= addrof ptr $MCC_GCTIB__LTest__A_3B, 8= addrof ptr $__classinforo__LTest__A_3B, 9= addrof ptr $classStateInitialized]
var $__cinf_LTest__B_3B <$__class_meta__> public
var $__fields_info__LTest__B_3B fstatic <[1] <$__field_info__>> public = [[1= 4, 2= 0, 3= 0xfdc0, 4= 0, 5= 232, 6= 352, 7= 172, 8= addrof ptr $__cinf_LTest__B_3B]]
var $__methods_info__LTest__B_3B fstatic <[1] <$__method_info__>> public = [[1= 0xfff6, 2= addrof ptr $__cinf_LTest__B_3B, 3= addroffunc ptr &LTest__B_3B_7C_3Cinit_3E_7C_28_29V, 4= 0x10000, 5= 188, 6= 216, 7= 172, 8= 0xad81, 9= 1, 10= 0]]
var $__superclasses__LTest__B_3B fstatic <[1] <$__superclass_meta__>> public = [[1= 0x4000000000000000]]
var $__classinforo__LTest__B_3B fstatic <$__class_meta_ro__> public = [1= 280, 2= addrof ptr $__fields_info__LTest__B_3B, 3= addrof ptr $__methods_info__LTest__B_3B, 4= addrof ptr $__superclasses__LTest__B_3B, 5= 1, 6= 1, 7= 0, 8= 1, 9= 0, 10= 32, 11= 172, 12= 0]
var $MCC_GCTIB__LTest__B_3B fstatic <* void> public
var $__cinf_LTest__B_3B <$__class_meta__> public = [1= 0x13bdf92a2d6fcf1e, 2= 0, 3= 0xffff, 4= 0, 5= 0, 6= addrof ptr $__vtb_LTest__B_3B, 7= addrof ptr $MCC_GCTIB__LTest__B_3B, 8= addrof ptr $__classinforo__LTest__B_3B, 9= addrof ptr $classStateInitialized]
var $__muid_classmetadata_bucket$$RCCycleTest_jar <[3] <* void>> public = [addrof ptr $__cinf_LRCCycleTest_3B, addrof ptr $__cinf_LTest__A_3B, addrof ptr $__cinf_LTest__B_3B]
func &MCC_Reflect_ThrowCastException nosideeffect () void
func &MCC_Reflect_Check_Casting_NoArray nosideeffect () void
func &MCC_Reflect_Check_Casting_Array nosideeffect () void
func &MCC_CheckThrowPendingException nosideeffect () void
func &MCC_PreNativeCall (var %caller ref) <* void>
func &MCC_PostNativeCall (var %env <* void>) void
func &MCC_DecodeReference nosideeffect (var %obj ref) ref
func &MCC_CallFastNative (var %func <* void>) <* void>
func &MCC_CallSlowNative0 (var %func <* void>) <* void>
func &MCC_CallSlowNative1 (var %func <* void>) <* void>
func &MCC_CallSlowNative2 (var %func <* void>) <* void>
func &MCC_CallSlowNative3 (var %func <* void>) <* void>
func &MCC_CallSlowNative4 (var %func <* void>) <* void>
func &MCC_CallSlowNative5 (var %func <* void>) <* void>
func &MCC_CallSlowNative6 (var %func <* void>) <* void>
func &MCC_CallSlowNative7 (var %func <* void>) <* void>
func &MCC_CallSlowNative8 (var %func <* void>) <* void>
func &MCC_CallFastNativeExt (var %func <* void>) <* void>
func &MCC_CallSlowNativeExt (var %func <* void>) <* void>
func &MCC_SetReliableUnwindContext nosideeffect () void
var $__reg_jni_func_tab$$RCCycleTest_jar <[0] <* void>>
var $__cinf_Ljava_2Flang_2FSystem_3B extern <$__class_meta__>
func &MCC_getFuncPtrFromItabSecondHash64 nosideeffect () ptr
var $__muid_func_def_tab$$RCCycleTest_jar fstatic <[5] <$MUIDFuncDefTabEntry>> = [[1= addroffunc ptr &LRCCycleTest_3B_7C_3Cinit_3E_7C_28_29V], [1= addroffunc ptr &LRCCycleTest_3B_7Cmain_7C_28ALjava_2Flang_2FString_3B_29V], [1= addroffunc ptr &LTest__A_3B_7C_3Cinit_3E_7C_28_29V], [1= addroffunc ptr &LTest__A_3B_7CgenCycle_7C_28_29V], [1= addroffunc ptr &LTest__B_3B_7C_3Cinit_3E_7C_28_29V]]
var $__muid_func_inf_tab$$RCCycleTest_jar fstatic <[5] <$MUIDFuncInfTabEntry>> = [[1= addroffunc ptr &LRCCycleTest_3B_7C_3Cinit_3E_7C_28_29V, 2= addroffunc ptr &LRCCycleTest_3B_7C_3Cinit_3E_7C_28_29V], [1= addroffunc ptr &LRCCycleTest_3B_7Cmain_7C_28ALjava_2Flang_2FString_3B_29V, 2= addroffunc ptr &LRCCycleTest_3B_7Cmain_7C_28ALjava_2Flang_2FString_3B_29V], [1= addroffunc ptr &LTest__A_3B_7C_3Cinit_3E_7C_28_29V, 2= addroffunc ptr &LTest__A_3B_7C_3Cinit_3E_7C_28_29V], [1= addroffunc ptr &LTest__A_3B_7CgenCycle_7C_28_29V, 2= addroffunc ptr &LTest__A_3B_7CgenCycle_7C_28_29V], [1= addroffunc ptr &LTest__B_3B_7C_3Cinit_3E_7C_28_29V, 2= addroffunc ptr &LTest__B_3B_7C_3Cinit_3E_7C_28_29V]]
var $__muid_func_def_muid_tab$$RCCycleTest_jar fstatic <[5] <$MUIDFuncDefMuidTabEntry>> = [[1= 0x53f926f4d5821784, 2= -849741229411409875], [1= -2057626659996252319, 2= -4244148982706990495], [1= -6522593282547121157, 2= -589866838482586098], [1= 0x2550930ad0ab5597, 2= -1656413819215847323], [1= 0x749dff24bd4960cd, 2= -1856617688624824905]]
var $__muid_func_muid_idx_tab$$RCCycleTest_jar fstatic <[5] u32> = [1, 4, 3, 0, 2]
var $__muid_data_def_tab$$RCCycleTest_jar fstatic <[3] <$MUIDDataDefTabEntry>> = [[1= addrof ptr $__cinf_LTest__B_3B], [1= addrof ptr $__cinf_LRCCycleTest_3B], [1= addrof ptr $__cinf_LTest__A_3B]]
var $__muid_data_def_muid_tab$$RCCycleTest_jar fstatic <[3] <$MUIDDataDefMuidTabEntry>> = [[1= -2673171647261450450, 2= -3219846449415409914], [1= -5627869524519406828, 2= -1366657250152391875], [1= 0x6aec35bbdfcf1c7f, 2= -867112758802146249]]
var $__muid_func_undef_tab$$RCCycleTest_jar fstatic <[12] <$MUIDUnifiedUndefTabEntry>> = [[1= 0], [1= 0], [1= 0], [1= 0], [1= 0], [1= 0], [1= 0], [1= 0], [1= 0], [1= 0], [1= 0], [1= 0]]
var $__muid_func_undef_muid_tab$$RCCycleTest_jar fstatic <[12] <$MUIDUnifiedUndefMuidTabEntry>> = [[1= 0x3e32352aee789835, 2= -3887705395317205813], [1= -4187412136968710015, 2= -3800091941095621250], [1= 0x6742c234127e0a27, 2= -3762262047879347071], [1= 0x7230554331c55d92, 2= -3676689525926909155], [1= -715372855679083712, 2= -2647497990906227723], [1= -7464356948810446352, 2= -2259485500590180091], [1= 0x783627f2afd1cbde, 2= -2046851302095768916], [1= -2701934576591406938, 2= -1693831364093527548], [1= 0x126777a7fe39e1fb, 2= -1314856249532362766], [1= 0xbf40578f3343f7a, 2= -1198421541845410999], [1= 0x7ca2bdf69a6c7c94, 2= -801329978528900548], [1= 0x477aafa4d7dd102b, 2= -442561182569419835]]
var $__muid_data_undef_tab$$RCCycleTest_jar fstatic <[3] <$MUIDUnifiedUndefTabEntry>> = [[1= addrof ptr $__cinf_Ljava_2Flang_2FObject_3B], [1= addrof ptr $Ljava_2Flang_2FSystem_3B_7Cout], [1= addrof ptr $__cinf_Ljava_2Flang_2FSystem_3B]]
var $__muid_data_undef_muid_tab$$RCCycleTest_jar fstatic <[3] <$MUIDUnifiedUndefMuidTabEntry>> = [[1= -567417612161374449, 2= -3298852447504547670], [1= 0x191283ac418c4bb9, 2= -1676204161023949463], [1= -5921653145571052587, 2= -171150348656858163]]
var $__muid_range_tab$$RCCycleTest_jar fstatic <[29] <$MUIDRangeTabEntry>> = [[1= -8525906493570096523, 2= -712220437341393984], [1= 0x6e52ae1238c0bb48, 2= -2364429413554491447], [1= 2, 2= 2], [1= 3, 2= 3], [1= 4, 2= 4], [1= 5, 2= 5], [1= 6, 2= 6], [1= 7, 2= 7], [1= 8, 2= 8], [1= 9, 2= 9], [1= 0, 2= 0], [1= 11, 2= 11], [1= 12, 2= 12], [1= 13, 2= 13], [1= 14, 2= 14], [1= 15, 2= 15], [1= addrof ptr $__muid_func_def_tab$$RCCycleTest_jar, 2= addrof ptr $__muid_func_def_tab$$RCCycleTest_jar], [1= 0, 2= 0], [1= addrof ptr $__muid_func_inf_tab$$RCCycleTest_jar, 2= addrof ptr $__muid_func_inf_tab$$RCCycleTest_jar], [1= addrof ptr $__muid_func_undef_tab$$RCCycleTest_jar, 2= addrof ptr $__muid_func_undef_tab$$RCCycleTest_jar], [1= addrof ptr $__muid_data_def_tab$$RCCycleTest_jar, 2= addrof ptr $__muid_data_def_tab$$RCCycleTest_jar], [1= 0, 2= 0], [1= addrof ptr $__muid_data_undef_tab$$RCCycleTest_jar, 2= addrof ptr $__muid_data_undef_tab$$RCCycleTest_jar], [1= addrof ptr $__muid_func_def_muid_tab$$RCCycleTest_jar, 2= addrof ptr $__muid_func_def_muid_tab$$RCCycleTest_jar], [1= addrof ptr $__muid_func_undef_muid_tab$$RCCycleTest_jar, 2= addrof ptr $__muid_func_undef_muid_tab$$RCCycleTest_jar], [1= addrof ptr $__muid_data_def_muid_tab$$RCCycleTest_jar, 2= addrof ptr $__muid_data_def_muid_tab$$RCCycleTest_jar], [1= addrof ptr $__muid_data_undef_muid_tab$$RCCycleTest_jar, 2= addrof ptr $__muid_data_undef_muid_tab$$RCCycleTest_jar], [1= addrof ptr $__muid_func_muid_idx_tab$$RCCycleTest_jar, 2= addrof ptr $__muid_func_muid_idx_tab$$RCCycleTest_jar], [1= 0, 2= 0]]
var $__reflection_strtab$$RCCycleTest_jar fstatic <[91] u8> = [0, 76, 82, 67, 67, 121, 99, 108, 101, 84, 101, 115, 116, 59, 0, 109, 97, 105, 110, 0, 40, 91, 76, 106, 97, 118, 97, 47, 108, 97, 110, 103, 47, 83, 116, 114, 105, 110, 103, 59, 41, 86, 0, 48, 33, 48, 0, 60, 105, 110, 105, 116, 62, 0, 40, 41, 86, 0, 76, 84, 101, 115, 116, 95, 65, 59, 0, 98, 98, 0, 76, 84, 101, 115, 116, 95, 66, 59, 0, 103, 101, 110, 67, 121, 99, 108, 101, 0, 97, 97, 0]
var $__compilerVersionNum$$RCCycleTest_jar <[0] <* void>> = [1, 0]
func &LRCCycleTest_3B_7C_3Cinit_3E_7C_28_29V public constructor (var %_this <* <$LRCCycleTest_3B>>) void {
  funcid 48153
  var %Reg1_R43696 <* <$LRCCycleTest_3B>> localrefvar
  var %Reg1_R57 <* <$Ljava_2Flang_2FObject_3B>> localrefvar
  var %__muid_symptr <* void>

  intrinsiccall MCCIncRef (dread ref %_this)
  intrinsiccall MCCDecRef (dread ref %Reg1_R43696)
  dassign %Reg1_R43696 0 (dread ref %_this)
  #INSTIDX : 0||0000:  aload_0
  #INSTIDX : 1||0001:  invokespecial
  regassign ptr %1 (dread ref %Reg1_R57)
  dassign %Reg1_R57 0 (retype ref <* <$Ljava_2Flang_2FObject_3B>> (dread ref %Reg1_R43696))
  intrinsiccall MCCIncRef (dread ref %Reg1_R57)
  intrinsiccall MCCDecRef (regread ptr %1)
  #Call function:Ljava_2Flang_2FObject_3B_7C_3Cinit_3E_7C_28_29V
  dassign %__muid_symptr 0 (iread ptr <* <$MUIDUnifiedUndefTabEntry>> 1 (array 0 ptr <* <[12] <$MUIDUnifiedUndefTabEntry>>> (addrof ptr $__muid_func_undef_tab$$RCCycleTest_jar, constval i64 9)))
  icallassigned (dread ptr %__muid_symptr, dread ref %Reg1_R57) {}
  #INSTIDX : 4||0004:  return
  intrinsiccall MPL_CLEANUP_LOCALREFVARS (dread ref %Reg1_R43696, dread ref %Reg1_R57)
  return ()
}
func &LRCCycleTest_3B_7Cmain_7C_28ALjava_2Flang_2FString_3B_29V public static (var %Reg3_R743 <* <[] <* <$Ljava_2Flang_2FString_3B>>>>) void {
  funcid 48154
  var %Reg0_R43699 <* <$LTest__A_3B>> localrefvar
  var %Reg2_R43699 <* <$LTest__A_3B>> localrefvar
  var %Reg0_R43698 <* <$LTest__B_3B>> localrefvar
  var %Reg0_R562 <* <$Ljava_2Fio_2FPrintStream_3B>> localrefvar
  var %Reg1_R43 <* <$Ljava_2Flang_2FString_3B>> localrefvar
  var %L_STR_161337 <* <$Ljava_2Flang_2FString_3B>>
  var %L_STR_161338 <* <$Ljava_2Flang_2FString_3B>>
  var %__muid_symptr <* void>

  #INSTIDX : 0||0000:  new
  regassign ptr %1 (dread ref %Reg0_R43699)
  dassign %Reg0_R43699 0 (gcmalloc ref <$LTest__A_3B>)
  intrinsiccall MCCDecRef (regread ptr %1)
  #INSTIDX : 3||0003:  dup
  #INSTIDX : 4||0004:  invokespecial
  #Call function:LTest__A_3B_7C_3Cinit_3E_7C_28_29V
  dassign %__muid_symptr 0 (iread ptr <* <$MUIDFuncDefTabEntry>> 1 (array 0 ptr <* <[5] <$MUIDFuncDefTabEntry>>> (addrof ptr $__muid_func_def_tab$$RCCycleTest_jar, constval i64 2)))
  icallassigned (dread ptr %__muid_symptr, dread ref %Reg0_R43699) {}
  #INSTIDX : 7||0007:  astore_1
  intrinsiccall MCCIncRef (dread ref %Reg0_R43699)
  intrinsiccall MCCDecRef (dread ref %Reg2_R43699)
  dassign %Reg2_R43699 0 (dread ref %Reg0_R43699)
  #INSTIDX : 8||0008:  aload_1
  #INSTIDX : 9||0009:  invokevirtual
  icallassigned (
    iread u64 <* u64> 0 (add ptr (
      iread ptr <* <$__class_meta__>> 6 (iread ref <* <$Ljava_2Flang_2FObject_3B>> 1 (dread ref %Reg2_R43699)),
      constval u32 88)),
    dread ref %Reg2_R43699) {}
  #INSTIDX : 12||000c:  aload_1
  #INSTIDX : 13||000d:  getfield
  intrinsiccallassigned MCCLoadRef (
    dread ref %Reg2_R43699,
    iaddrof ref <* <$LTest__A_3B>> 4 (dread ref %Reg2_R43699)) { regassign ptr %2}
  intrinsiccall MCCDecRef (dread ref %Reg0_R43698)
  dassign %Reg0_R43698 0 (regread ptr %2)
  #INSTIDX : 16||0010:  ifnonnull
  brtrue @label0 (ne i32 ref (dread ref %Reg0_R43698, constval ref 0))
  #INSTIDX : 19||0013:  getstatic
  intrinsiccall MPL_CLINIT_CHECK (addrof ptr $__cinf_Ljava_2Flang_2FSystem_3B)
  regassign ptr %3 (dread ref %Reg0_R562)
  #Read from: Ljava_2Flang_2FSystem_3B_7Cout
  dassign %Reg0_R562 0 (iread ref <* <* <$Ljava_2Fio_2FPrintStream_3B>>> 0 (iread ptr <* <$MUIDUnifiedUndefTabEntry>> 1 (array 0 ptr <* <[3] <$MUIDUnifiedUndefTabEntry>>> (addrof ptr $__muid_data_undef_tab$$RCCycleTest_jar, constval i64 1))))
  intrinsiccall MCCIncRef (dread ref %Reg0_R562)
  intrinsiccall MCCDecRef (regread ptr %3)
  #INSTIDX : 22||0016:  ldc
  callassigned &MCC_GetOrInsertLiteral (addrof ptr $_C_STR_d4b6a4788952b3ae6554bbc29994ed26) { dassign %L_STR_161337 0 }
  intrinsiccall MCCIncRef (dread ptr %L_STR_161337)
  intrinsiccall MCCDecRef (dread ref %Reg1_R43)
  dassign %Reg1_R43 0 (dread ptr %L_STR_161337)
  #INSTIDX : 24||0018:  invokevirtual
  icallassigned (
    iread u64 <* u64> 0 (add ptr (
      iread ptr <* <$__class_meta__>> 6 (iread ref <* <$Ljava_2Flang_2FObject_3B>> 1 (dread ref %Reg0_R562)),
      constval u32 312)),
    dread ref %Reg0_R562,
    dread ref %Reg1_R43) {}
  #INSTIDX : 27||001b:  goto
  goto @label1
@label0   #INSTIDX : 30||001e:  getstatic
  intrinsiccall MPL_CLINIT_CHECK (addrof ptr $__cinf_Ljava_2Flang_2FSystem_3B)
  regassign ptr %4 (dread ref %Reg0_R562)
  #Read from: Ljava_2Flang_2FSystem_3B_7Cout
  dassign %Reg0_R562 0 (iread ref <* <* <$Ljava_2Fio_2FPrintStream_3B>>> 0 (iread ptr <* <$MUIDUnifiedUndefTabEntry>> 1 (array 0 ptr <* <[3] <$MUIDUnifiedUndefTabEntry>>> (addrof ptr $__muid_data_undef_tab$$RCCycleTest_jar, constval i64 1))))
  intrinsiccall MCCIncRef (dread ref %Reg0_R562)
  intrinsiccall MCCDecRef (regread ptr %4)
  #INSTIDX : 33||0021:  ldc
  callassigned &MCC_GetOrInsertLiteral (addrof ptr $_C_STR_a360f9874cf7de12fce9dffe2c4c052a) { dassign %L_STR_161338 0 }
  intrinsiccall MCCIncRef (dread ptr %L_STR_161338)
  intrinsiccall MCCDecRef (dread ref %Reg1_R43)
  dassign %Reg1_R43 0 (dread ptr %L_STR_161338)
  #INSTIDX : 35||0023:  invokevirtual
  icallassigned (
    iread u64 <* u64> 0 (add ptr (
      iread ptr <* <$__class_meta__>> 6 (iread ref <* <$Ljava_2Flang_2FObject_3B>> 1 (dread ref %Reg0_R562)),
      constval u32 312)),
    dread ref %Reg0_R562,
    dread ref %Reg1_R43) {}
@label1   #INSTIDX : 38||0026:  return
  intrinsiccall MPL_CLEANUP_LOCALREFVARS (dread ref %Reg0_R43699, dread ref %Reg2_R43699, dread ref %Reg0_R43698, dread ref %Reg0_R562, dread ref %Reg1_R43)
  return ()
}
func &LTest__A_3B_7C_3Cinit_3E_7C_28_29V constructor (var %_this <* <$LTest__A_3B>>) void {
  funcid 48155
  var %Reg1_R43699 <* <$LTest__A_3B>> localrefvar
  var %Reg1_R57 <* <$Ljava_2Flang_2FObject_3B>> localrefvar
  var %__muid_symptr <* void>

  intrinsiccall MCCIncRef (dread ref %_this)
  intrinsiccall MCCDecRef (dread ref %Reg1_R43699)
  dassign %Reg1_R43699 0 (dread ref %_this)
  #INSTIDX : 0||0000:  aload_0
  #INSTIDX : 1||0001:  invokespecial
  regassign ptr %1 (dread ref %Reg1_R57)
  dassign %Reg1_R57 0 (retype ref <* <$Ljava_2Flang_2FObject_3B>> (dread ref %Reg1_R43699))
  intrinsiccall MCCIncRef (dread ref %Reg1_R57)
  intrinsiccall MCCDecRef (regread ptr %1)
  #Call function:Ljava_2Flang_2FObject_3B_7C_3Cinit_3E_7C_28_29V
  dassign %__muid_symptr 0 (iread ptr <* <$MUIDUnifiedUndefTabEntry>> 1 (array 0 ptr <* <[12] <$MUIDUnifiedUndefTabEntry>>> (addrof ptr $__muid_func_undef_tab$$RCCycleTest_jar, constval i64 9)))
  icallassigned (dread ptr %__muid_symptr, dread ref %Reg1_R57) {}
  #INSTIDX : 4||0004:  return
  intrinsiccall MPL_CLEANUP_LOCALREFVARS (dread ref %Reg1_R43699, dread ref %Reg1_R57)
  return ()
}
func &LTest__A_3B_7CgenCycle_7C_28_29V virtual (var %_this <* <$LTest__A_3B>>) void {
  funcid 48156
  var %Reg3_R43699 <* <$LTest__A_3B>> localrefvar
  var %Reg0_R43698 <* <$LTest__B_3B>> localrefvar
  var %Reg2_R43698 <* <$LTest__B_3B>> localrefvar
  var %__muid_symptr <* void>

  intrinsiccall MCCIncRef (dread ref %_this)
  intrinsiccall MCCDecRef (dread ref %Reg3_R43699)
  dassign %Reg3_R43699 0 (dread ref %_this)
  #INSTIDX : 0||0000:  new
  regassign ptr %1 (dread ref %Reg0_R43698)
  dassign %Reg0_R43698 0 (gcmalloc ref <$LTest__B_3B>)
  intrinsiccall MCCDecRef (regread ptr %1)
  #INSTIDX : 3||0003:  dup
  #INSTIDX : 4||0004:  invokespecial
  #Call function:LTest__B_3B_7C_3Cinit_3E_7C_28_29V
  dassign %__muid_symptr 0 (iread ptr <* <$MUIDFuncDefTabEntry>> 1 (array 0 ptr <* <[5] <$MUIDFuncDefTabEntry>>> (addrof ptr $__muid_func_def_tab$$RCCycleTest_jar, constval i64 4)))
  icallassigned (dread ptr %__muid_symptr, dread ref %Reg0_R43698) {}
  #INSTIDX : 7||0007:  astore_1
  intrinsiccall MCCIncRef (dread ref %Reg0_R43698)
  intrinsiccall MCCDecRef (dread ref %Reg2_R43698)
  dassign %Reg2_R43698 0 (dread ref %Reg0_R43698)
  #INSTIDX : 8||0008:  aload_0
  #INSTIDX : 9||0009:  aload_1
  #INSTIDX : 10||000a:  putfield
  intrinsiccall MCCWrite (
    dread ref %Reg3_R43699,
    iaddrof ref <* <$LTest__A_3B>> 4 (dread ref %Reg3_R43699),
    dread ref %Reg2_R43698)
  #INSTIDX : 13||000d:  aload_1
  #INSTIDX : 14||000e:  aload_0
  #INSTIDX : 15||000f:  putfield
  intrinsiccall MCCWrite (
    dread ref %Reg2_R43698,
    iaddrof ref <* <$LTest__B_3B>> 4 (dread ref %Reg2_R43698),
    dread ref %Reg3_R43699)
  #INSTIDX : 18||0012:  return
  intrinsiccall MPL_CLEANUP_LOCALREFVARS (dread ref %Reg3_R43699, dread ref %Reg0_R43698, dread ref %Reg2_R43698)
  return ()
}
func &LTest__B_3B_7C_3Cinit_3E_7C_28_29V constructor (var %_this <* <$LTest__B_3B>>) void {
  funcid 48157
  var %Reg1_R43698 <* <$LTest__B_3B>> localrefvar
  var %Reg1_R57 <* <$Ljava_2Flang_2FObject_3B>> localrefvar
  var %__muid_symptr <* void>

  intrinsiccall MCCIncRef (dread ref %_this)
  intrinsiccall MCCDecRef (dread ref %Reg1_R43698)
  dassign %Reg1_R43698 0 (dread ref %_this)
  #INSTIDX : 0||0000:  aload_0
  #INSTIDX : 1||0001:  invokespecial
  regassign ptr %1 (dread ref %Reg1_R57)
  dassign %Reg1_R57 0 (retype ref <* <$Ljava_2Flang_2FObject_3B>> (dread ref %Reg1_R43698))
  intrinsiccall MCCIncRef (dread ref %Reg1_R57)
  intrinsiccall MCCDecRef (regread ptr %1)
  #Call function:Ljava_2Flang_2FObject_3B_7C_3Cinit_3E_7C_28_29V
  dassign %__muid_symptr 0 (iread ptr <* <$MUIDUnifiedUndefTabEntry>> 1 (array 0 ptr <* <[12] <$MUIDUnifiedUndefTabEntry>>> (addrof ptr $__muid_func_undef_tab$$RCCycleTest_jar, constval i64 9)))
  icallassigned (dread ptr %__muid_symptr, dread ref %Reg1_R57) {}
  #INSTIDX : 4||0004:  return
  intrinsiccall MPL_CLEANUP_LOCALREFVARS (dread ref %Reg1_R43698, dread ref %Reg1_R57)
  return ()
}


text:
(*
 Copyright (C) 2018 Susi Lehtola

 This Source Code Form is subject to the terms of the Mozilla Public
 License, v. 2.0. If a copy of the MPL was not distributed with this
 file, You can obtain one at http://mozilla.org/MPL/2.0/.
*)

(* type: gga_exc *)

(* x has a different definition in the Chachiyo paper  *)
chachiyo_x := x -> 2/9 * (Pi/3)**(1/3) * (2**(-1/3) * x):

(* equation 1 *)
chachiyo_f0 := x -> (3*x^2 + Pi^2*log(x+1)) / ((3*x + Pi^2)*log(x+1)):

chachiyo_f := x -> chachiyo_f0(chachiyo_x(x)):

f := (rs, z, xt, xs0, xs1) -> gga_exchange(chachiyo_f, rs, z, xs0, xs1):

ext: mpl


text:
#int fact (int n) {
#  if(n != 1)
#    return foo(n - 1);
#  else return 1;
#}

func $fact (
 var %n i32, var %p1 a64, var %p2 a64) i32 {
 if (ne i32 i32 (dread i32 %n, constval i32 1)) {
    intrinsiccall JSOP_ADD (
      dread a64 %p1, dread a64 %p2, dread i32 %n)
    return (regread i32 %%retval)}
   else {
     return (constval i32 1) } }

max_stars_repo_name: MapleSystem/OpenArkCompiler

max_issues_repo_name: MapleSystem/OpenArkCompiler

max_forks_repo_name: MapleSystem/OpenArkCompiler


text:
    {   Plan.mpl    }

    {   Aggregate production planning for 6 months   }

    {   COPYRIGHT (C) 1987-1989,  Maximal Software   }

TITLE
    Plan_for_x_years ;

DATA
    NrOfYears = 5?;
    Dec = 12;

INDEX
    i = 1..NrOfYears ;
    j = 1..Dec ;

DATA
    Demand[j] := 1000 (90,81,102,100,140,141,180,170,140,110,90,90);

DECISION
    Invt[i,j]
    Prod[i,j]
    Work[i,j]
    OvrTm[i,j]
    Hired[i,j] -> Hir
    Fired[i,j] -> Fir

MODEL

     MIN   ProductionCost -> ProdCost  =

           SUM(i,j:     8.80 Invt
                    +  64.33 Prod
                    +  17600 Work
                    +    160 OvrTm
                    +   8800 Hired
                    +  52800 Fired) ;

SUBJECT TO

  ProdInvt[i=1,j=1] -> PrIn :  Invt  =  20000         + Prod - Demand ;
  ProdInvt[i>1,j=1] -> PrIn :  Invt  =  Invt[i-1,Dec] + Prod - Demand ;
  ProdInvt[i,j>1]   -> PrIn :  Invt  =  Invt[j-1]     + Prod - Demand ;

  HiredFired[i=1,j=1] -> HiFi :  Work - 45             =  Hired - Fired ;
  HiredFired[i>1,j=1] -> HiFi :  Work - Work[i-1,Dec]  =  Hired - Fired ;
  HiredFired[i,j>1]   -> HiFi :  Work - Work[j-1]      =  Hired - Fired ;

  ProdHours[i,j] -> PrHo  : 170 Work + OvrTm  >  0.07 Prod ;
  OvrTmWork[i,j]  -> OvWo  : OvrTm  <  50% * 170 Work ;

BOUNDS

    MaxInvt   :  Invt[j<Dec]  <  90000 ;
    CloseInvt :  Invt[j=Dec]  =  20000 ;

    MaxWork   :  Work  <  60 ;
    MaxProd   :  Prod  <  170000 ;
END

max_stars_repo_name: BjarniMax/mpl-models

max_issues_repo_name: BjarniMax/mpl-models

max_forks_repo_name: BjarniMax/mpl-models


text:
(*
 Copyright (C) 2017 M.A.L. Marques
               2019 Susi Lehtola

 This Source Code Form is subject to the terms of the Mozilla Public
 License, v. 2.0. If a copy of the MPL was not distributed with this
 file, You can obtain one at http://mozilla.org/MPL/2.0/.
*)

(* type: gga_exc *)
(* replace: "my_dilog\(" -> "xc_dilogarithm(" *)

`diff/my_dilog` := proc(g, x) -log(1 - g)/g * diff(g, x) end proc:

gg99_a := 3^(1/4)/(2*sqrt(2)*Pi^(3/2)):
gg99_b := 4*sqrt(3)*Pi^3:

(* Equation 22 in the paper, i.e.
   the solution of x = 2*Pi*sinh(r)/(3*cosh(r))^(1/3) *)

gg99_r_branch1 := x -> arcsinh( (gg99_a * x * sqrt(x^2 + (gg99_b + sqrt(gg99_b^2 - x^6))^(2/3))) / (gg99_b + sqrt(gg99_b^2 - x^6))^(1/6) ):

(* The second branch is from Andrew Gilbert via email *)

gg99_r_branch2 := x-> arcsinh(sqrt(x^3*(3/gg99_b)*cos(arctan(sqrt(1/(gg99_b^2)*x^6-1))/3))):

(* Glue the pieces together. The min and max are required
   to avoid float exceptions *)
gg99_r := x -> my_piecewise3(x < gg99_b^(1/3),
  gg99_r_branch1(m_min(x, gg99_b^(1/3) - 1e-10)),
  gg99_r_branch2(m_max(x, gg99_b^(1/3) + 1e-10))
  ):

(* Equation 21 *)
gg99_f0 := r -> (Pi^2 - 12*r*log(1 + exp(-2*r)) +
12*my_dilog(-exp(-2*r))) / (2*3^(1/3)*Pi*r*sech(r)^(2/3)) / X_FACTOR_C:

(* Assemble the function *)
gg99_f := x -> gg99_f0(gg99_r(x)):

f := (rs, zeta, xt, xs0, xs1) ->
  gga_exchange(gg99_f, rs, zeta, xs0, xs1):


text:
$ifndef _GEN_SOL_
$define _GEN_SOL_

$include "Utils.mpl"

# 寻找一般解,考虑存在常数约束和多个单变量非零约束
# 若只存在一个极大元，则直接返回下标
# 否则分组返回极大元，相等的极大元放在一组中
findGenSolInd:=proc(icons::list(set),n::posint)
    return findGenVec(con2vec~(icons,n));
end proc:

# 将约束条件转化为向量
con2vec:=proc(c::set,n::posint)
    local v,t;
    v:=Array(1..n,2);
    for t in c do
        if op(0,t)=`=` then
            v[op([1,1],indets(t,name))]:=0;
        elif op(0,t)=`<>` then
            v[op([1,1],indets(t,name))]:=1;
        end if;
    end do;
    return convert(v,list);
end proc:

# 寻找一般解
# 等价于寻找偏序集的极大元
# 返回下标的list
findGenVec:=proc(v::list)
    local mInd,i,j,n,isAdd,ind,t,gen,mgen,maxind;
    n:=numelems(v);
    mInd:=table();
    mInd[1]:=1;
    for i from 2 to n do
        isAdd:=true;
        for j in indices(mInd,nolist) do
            if v[i] &> v[j] then
                unassign(evaln(mInd[j]));
            elif v[j] &> v[i]  then
                isAdd:=false;
                break;
            end if;
        end do;
        if isAdd then
            mInd[i]:=1;
        end if;
    end do;
    ind:=[indices(mInd,nolist)];
    if numelems(ind)=1 then
        return ind[1];
    else
        # 分组返回，相等的元素放在一组
        t:=table();
        for i in ind do
            tappend(t,v[i],i);
        end do;
        ind:=[entries(t,nolist)];
        # 选择一般性最大的分组
        gen:=map(x->add(v[x[1]]),ind);
        mgen:=max(gen);
        maxind:=find(gen,mgen);
        ind:=ind[maxind];
        return ind;
    end if;
end proc:

`&>=`:=proc(x::list,y::list)
    local n:=numelems(x);
    return andmap(k->evalb(x[k]>=y[k]),[seq(1..n)]);
end proc:

`&>`:=proc(x::list,y::list)
    local n:=numelems(x);
    return x &>= y and ormap(k->evalb(x[k]>y[k]),[seq(1..n)]);
end proc:

$endif


text:
######################################################################
# partitions(A) is the set of all partitions of A
# (Section sec-partitions)

# Partitions are represented as nonempty sets of nonempty subsets, so
# {{1,2},{3},{4,5}} is an element of partitions({1,2,3,4,5}).

`is_element/partitions` := (A::set) -> proc(pi)
 local n,T,U;
 global reason;

 if not type(pi,set) then
  reason := [convert(procname,string),"pi is not a set",pi];
  return false;
 fi;

 n := nops(pi);
 if n = 0 then
  reason := [convert(procname,string),"pi is empty"];
  return false;
 fi;

 T := {};

 for U in pi do
  if U = {} then
   reason := [convert(procname,string),"pi contaisn the empty set",pi];
   return false;
  fi;

  if U minus A <> {} then
   reason := [convert(procname,string),"U is not contained in A",U,A];
   return false;
  fi;

  if T intersect U <> {} then
   reason := [convert(procname,string),"sets in pi are not disjoint",U,T];
   return false;
  fi;

  T := T union U;
 od;

 if A minus T <> {} then
  reason := [convert(procname,string),"sets in pi do not cover A",T,A];
  return false;
 fi;

 return true;
end;

######################################################################

`is_equal/partitions` := (A::set) -> (pi,rho) -> evalb(pi = rho):

######################################################################

`is_leq/partitions` := (A::set) -> proc(pi,rho)
 local T,U,ok;
 global reason;

 for T in rho do
  ok := false;
  for U in pi do
   if T minus U = {} then
    ok := true;
    break;
   fi;
  od;

  if not ok then
   reason := [convert(procname,string),"T is not contained in any block of pi",T,pi];
   return false;
  fi;
 od;

 return true;
end;

######################################################################

`random_element/partitions` := (A::set) -> proc()
 local n,m,r,t,i,a,pi;

 if A = {} then return FAIL; fi;

 n := nops(A);
 m := rand(1..n)();
 r := rand(1..m);
 t := table();
 for i from 1 to m do t[i] := NULL; od;
 for a in A do
  i := r();
  t[i] := t[i],a;
 od:

 pi := NULL;
 for i from 1 to m do
  if t[i] <> NULL then pi := pi,{t[i]}; fi;
 od;

 pi := {pi};

 return pi;
end;

######################################################################

`list_elements/partitions` := proc(A::set)
 local n,a,A0,B0,B,C,S,P,Q,rho,pi;

 n := nops(A);
 if n = 0 then return {}; fi;
 if n = 1 then return {{A}}; fi;

 A0 := sort([op(A)]);
 a := A0[1];
 A0 := {op(2..-1,A0)};
 P := {A};
 S := `list_elements/subsets`(A0);

 for B0 in S do
  B := {a,op(B0)};
  Q := `list_elements/partitions`(A minus B);
  for rho in Q do
   pi := {B,op(rho)};
   P := P,pi;
  od;
 od;

 P := [P];
 return P;
end;

######################################################################

`count_elements/partitions` := (A::set) -> combinat[bell](nops(A));

`top/partitions` := (A::set) -> `if`(A = {},FAIL,map(a -> {a},A));

`bot/partitions` := (A::set) -> `if`(A = {},FAIL,{A});


text:
(*
 Copyright (C) 2017 M.A.L. Marques
               2020 Susi Lehtola

 This Source Code Form is subject to the terms of the Mozilla Public
 License, v. 2.0. If a copy of the MPL was not distributed with this
 file, You can obtain one at http://mozilla.org/MPL/2.0/.
*)

(* type: gga_exc *)

pbepow_kappa := KAPPA_PBE:
pbepow_mu    := 0.2195149727645171:

(* pbepow_m has to be written in float format otherwise the code is converted to C
   in a numerically unstable way *)
pbepow_m     := 100.0:

(* The geometrical series 1 + pbepow_kappa/pbepow_m * add(pbepow_yy(X2S*s)^i, i=1..pbepow_m)
   is summed analytically and simplified *)
pbepow_f0 := s -> 1 - pbepow_mu*s^2 * ( ((pbepow_m * pbepow_mu * s^2) / (pbepow_m * pbepow_mu * s^2 + pbepow_kappa))^pbepow_m - 1):

pbepow_f  := x -> pbepow_f0(X2S*x):

f := (rs, z, xt, xs0, xs1) -> gga_exchange(pbepow_f, rs, z, xs0, xs1):


text:
SemiAlgebraic := module()
    uses Domain, SolveTools;
    export SimplName  := "SemiAlgebraic";
    export SimplOrder := 6+(1/2);

    export ModuleApply := proc(dbnds :: DomBound, dshape :: DomShape, $)
        local vs;
        vs := ListTools[Reverse](Domain:-Bound:-varsOf(dbnds));
        vs := map(v->if v::list then op(1,v) else v end if, vs);
        do_sh( dshape , dbnds, vs );
    end proc;

    local extra_sol := (c -> c::`=` and (lhs(c)::DomBoundVar and evalb(lhs(c)=rhs(c))));
    local postproc := proc(sol, $)
      foldl((x,k)->subsindets(x,op(k)), sol
           ,[specfunc('piecewise')
            ,x->DSplit(Partition:-PWToPartition(x, 'do_solve'))]
           ,[list(list),DSum@op]
           ,[list(relation),x->DConstrain(remove(extra_sol,x)[])]);
    end proc;

    local do_sh := proc( sh :: DomShape, ctx :: DomBound, vs, $)
        local sol;
        if sh :: DomConstrain then
            sol := do_Constrain(sh, ctx, vs);
            if sol :: DomShape then sol
            else postproc(sol) end if;
        elif sh :: DomSplit then
            DSplit(Partition:-Amap(
              [ (c,_)->c
              , (p,c)-> p->do_sh(p, applyop(cs->{c,op(cs)},2,ctx), vs)
              , c->c ], op(1, sh) ));
        elif sh :: DomSum then
            map(s->do_sh(s, ctx, vs), sh);
        else
            DNoSol(sprintf("Don't know how to solve DOMAIN(%a, %a)", ctx, sh));
        end if;
    end proc;

    local do_Constrain := proc( sh :: DomConstrain , ctx, vs_, $ )
        local ts, vs := vs_, cs;
        cs := { op( Domain:-Bound:-toConstraints(ctx,'no_infinity') ), op(sh) } ;
        if nops(ctx) >= 2 then cs := cs union op(2,ctx) end if;
        ts, cs := selectremove(type, cs, `::`);
        ts := select(t->op(2,t)<>real, ts);
        if ts <> {} and has(cs,ts) then
          WARNING("Variables of unrecognized types being passed to SemiAlgebraic: %1",ts);
        end if;
        try SemiAlgebraic(cs, vs);
        catch: DNoSol(StringTools[FormatMessage](lastexception[2..-1]));
        end try;
    end proc;
end module; #SemiAlgebraic

max_stars_repo_name: zaxtax/hakaru

max_issues_repo_name: zaxtax/hakaru

max_forks_repo_name: zaxtax/hakaru


text:
# Base Parameter Energy Regressor for Robot based on MDH frames
# Einleitung
# Erstellung einer parameterlinearen und -minimalen Regressorform
#
# Dateiname:
# robot -> Berechnung für allgemeinen Roboter
# chain -> Berechnung für eine serielle Struktur (nicht: Baumstruktur)
# fixb -> fixed base. Kein Floating base Modell. Dort ist diese Form der Minimalparameterform nicht möglich.
# energy -> Berechnung bezogen auf Energie
# regressor -> Regressorform (parameterlinear)
# linearsolve -> Nicht geometrischer Ansatz nach Khalil, sondern symbolischer Ansatz mit LinearSolve
# Authors
# Jonas Diekmeyer (Studienarbeit bei Elias Knöchelmann)
# Moritz Schappler, moritz.schappler@imes.uni-hannover.de, 2020-06
#
# (C) Institut fuer Mechatronische Systeme, Leibniz Universitaet Hannover
# Quellen
# [Diekmeyer2018_S678] Identifikation der inversen Dynamik eines seriellen Roboters im geschlossenen Regelkreis, Studienarbeit, imes, LUH, 2018
# Initialisierung
interface(warnlevel=0): # Unterdrücke die folgende Warnung.
restart: # Gibt eine Warnung, wenn über Terminal-Maple mit read gestartet wird.
interface(warnlevel=3):
interface(rtablesize=100): # Zur Anzeige von größeren Vektoren
;
with(LinearAlgebra):
with(ArrayTools):
with(codegen):
with(CodeGeneration):
with(StringTools):
codegen_act := true:
codegen_opt := 2:
read "../helper/proc_convert_s_t":
read "../helper/proc_convert_t_s":
read "../helper/proc_MatlabExport":
read "../helper/proc_simplify2":
read "../robot_codegen_definitions/robot_env":
read sprintf("../codeexport/%s/tmp/tree_floatb_definitions", robot_name, base_method_name):
read sprintf("../codeexport/%s/tmp/kinematic_constraints_maple_inert.m", robot_name):
kin_constraints_exist := kin_constraints_exist: # nur zum Abschätzen der Komplexität
;
# Prüfe, ob die symbolische Berechnung der Parameterminimierung berechnet werden sollte
abort_this_worksheet := false:
# Bestimme, ob es eine Baumstruktur ist. Wenn ja, funktioniert der andere Algorithmus nicht und dieser wird genommen.
tree:=false:
for i from 1 to NJ do
  if not v(i) = i-1 then
    tree := true: break:
  end if:
end do:

if not (assigned(user_CoM) or assigned(user_M) or assigned(user_inertia) \
  or kin_constraints_exist or tree) then
  # Es gibt keinen Sonderfall, diese Berechnung der Minimalparameter ist nicht notwendig.
  printf("Keine analytische Berechnung der Minimalparameter notwendig.\n"):
  abort_this_worksheet := true:
end if:

if assigned(dynpar_minimization_linsolve) and not dynpar_minimization_linsolve then
  printf("Analytische Berechnung der Minimalparameter manuell deaktiviert.\n"):
  abort_this_worksheet := true:
end if:

if assigned(dynpar_minimization_linsolve) and dynpar_minimization_linsolve then
  printf("Analytische Berechnung der Minimalparameter durch Option dynpar_minimization_linsolve verlangt.\n"):
  abort_this_worksheet := false:
end if:

if abort_this_worksheet then
  quit: # Funktioniert in GUI nicht richtig...
  robot_name := "": # ...Daher auch Löschung des Roboternamens.
end if:
printf("%s. Generiere Minimalparameterregressor der Energie für %s (symbolischer Ansatz)\n", \
  FormatTime("%Y-%m-%d %H:%M:%S"), robot_name, codegen_dynpar):
# Term-Vereinfachungen einstellen
if not assigned(simplify_options) or simplify_options(7)=-1 then # Standard-Einstellungen:
  if not kin_constraints_exist then # normale serielle Ketten und Baumstrukturen
    use_simplify := 0: # Standardmäßig aus
  else # mit kinematischen Zwangsbedingungen
    use_simplify := 1: # standardmäßig simplify-Befehle anwenden
  end if:
else # Benutzer-Einstellungen:
  use_simplify := simplify_options(7): # siebter Eintrag ist für Energie-Regressor
end if:

# Ergebnisse der Energie laden (aus robot_chain_fixb_rotmat_energy_regressor.mw)
read sprintf("../codeexport/%s/tmp/energy_kinetic_fixb_regressor_maple.m", robot_name):
t_ges := t_ges:
read sprintf("../codeexport/%s/tmp/energy_potential_fixb_regressor_maple.m", robot_name):
u_ges := u_ges:
# Parameterminimierung
# Siehe auch [Diekmeyer2018_S678] S. 19f. und 32ff.
# Funktionen definieren

splitSummands := proc(input_expr)
# Gibt eine Liste aller Summanden zurücachezurück
  local expr, sub_expr:

  expr := expand(input_expr):

  if type(expr,'`+`') then
    sub_expr := op(expr):
    return sub_expr:
  else
    return expr:
  end if:

end proc:

splitFactors := proc(input_expr)
# Gibt eine Liste aller Faktoren zurück
  local sub_expr:

  if type(input_expr,'`*`') then
    sub_expr := op(input_expr):
    return sub_expr:
  else
    return input_expr:
  end if:

end proc:

removeEmptyRows := proc(h,G)
  local mh, mG, i, j, flag:

  flag := Vector[column](Size(G,1)):
  for i from 1 to Size(G,1) do
    flag(i) := 0:
    for j from 1 to Size(G,2) do
      if G(i,j) <> 0 then
        flag(i) := 1:
      end if:
    end do:
  end do:

  mh := Vector[column](Norm(flag,1)):
  mG := Matrix(Norm(flag,1),Size(G,2)):

  j := 1:
  for i from 1 to Size(G,1) do
    if flag(i) = 1 then
      mG(j,..) := G(i,..):
      mh(j) := h(i):
      j := j + 1:
    end if:
  end do:

  return mh, mG:
end proc:

permuteEmptyRows := proc(G)
  local Pb, Pd, i, j, k, flag:

  flag := Vector[column](Size(G,1)):
  for i from 1 to Size(G,1) do
    flag(i) := 0:
    for j from 1 to Size(G,2) do
      if G(i,j) <> 0 then
        flag(i) := 1:
      end if:
    end do:
  end do:

  Pb := Matrix(Norm(flag,1),Size(G,1)):
  Pd := Matrix(Size(G,1)-Norm(flag,1),Size(G,1)):

  j := 0:
  k := 0:
  for i from 1 to Size(G,1) do
    if flag(i) = 1 then
      j := j + 1:
      Pb(j,i) := 1:
    elif flag(i) = 0 then
      k := k + 1:
      Pd(k,i) := 1:
    end if:
  end do:

  return Pb, Pd:
end proc:

hasElement := proc(v, expr)
# Gibt den Index des Vektorelements von "v" wieder, das dem Ausdruck "expr" entspricht.
# Ist der Ausdruck nicht enthalten wird der höchste Index um eins erhöht wiedergegeben.
  local i:

  for i from 1 to NumElems(v) do
    if verify(v(i),expr) then
      return i:
    end if:
  end do:

  return NumElems(v)+1:
end proc:

tmp_tt1 := time(): tmp_timelastmessage := tmp_tt1:
# Kinematikparameter, die als Argument von Sinus und Cosinus stehen ersetzen.
# Ansonsten werden sie durch die Umwandlung in die Exponentialform mit ersetzt. Das ist aber nicht hilfreich.
read sprintf("../codeexport/%s/tmp/parameter_kin", robot_name):
pkin := pkin:

lagrange := Vector[column](t_ges[1,..] - u_ges[1,..]):
(*
# TODO: Das scheint die Rechenzeit von LinearSolve in manchen Fällen stark zu verlängern
for j from 1 to RowDimension(pkin) do
  lagrange := subs(parse(sprintf("sin(%s)", pkin(j,1)))=parse(sprintf("sin_%s", pkin(j,1))), lagrange):
  lagrange := subs(parse(sprintf("cos(%s)", pkin(j,1)))=parse(sprintf("cos_%s", pkin(j,1))), lagrange):
end do:
*)

# Parameter gruppieren. Siehe [Diekmeyer2018_S678] Gl. 3.27
# Erstelle eine Matrix V, die die einzelnen Terme des Energie-Regressors gruppiert
# Komponentenmatrix für jeden Eintrag des Lagrange-Regressors.
# Initialisierung mit voller Dimension (notwendig für geschlossene kinematische Ketten;
# in dem Fall sind die Einträge der letzten Parameter Null/konstant. Bei rein seriellen Ketten nicht der Fall).
V := Matrix(NumElems(lagrange),1):
w := Vector[column](): # Hilfsvariable
for i from 1 to NumElems(lagrange) do # Parameter-Einträge durchgehen
  if i=1 or time()-tmp_timelastmessage > 10.0 then # min. alle 10s Status ausgeben.
    printf("%s. Erstelle Komponentenmatrix für Lagrange-Regressor-Eintrag %d/%d\n", FormatTime("%Y-%m-%d %H:%M:%S"), i, NumElems(lagrange)):
    tmp_timelastmessage := time():
  end if:
  lagrange_term := simplify(expand(lagrange(i))):
  # Umrechnung der Trigonometrische Funktion in Exponent-Darstellung.
  # Funktioniert nicht ganz wie erwartet, falls Parameter in sin/cos stehen. Dann treten Imaginärteile auf.
  # Diese werden weiter unten wieder entfernt
  # Durch diesen Befehl wird die Spalten-Dimension der Matrix stark reduziert (ca. auf ein Drittel)
  lagrange_term := convert(lagrange_term,exp):
  lagrange_term := simplify(expand(lagrange_term),size):
  lagrange_summand := [splitSummands(lagrange_term)]:
# print(lagrange_summand):
  for j from 1 to nops(lagrange_summand) do
#   printf("j=%d\n", j):
    # Konstante Ausdrücke werden nicht berücksichtigt!
    if not has(lagrange_summand[j],t) then
      next:
    end if:
    lagrange_factor := [splitFactors(lagrange_summand[j])]:
    lagrange_timeVar := 1:
    lagrange_constant := 1:
    for k from 1 to nops(lagrange_factor) do
#     printf("k=%d\n", k):
      if has(lagrange_factor[k],t) then
        lagrange_timeVar := lagrange_timeVar * lagrange_factor[k]:
      else
        lagrange_constant := lagrange_constant * lagrange_factor[k]:
      end if
    end do:
    lagrange_timeVar := simplify(expand(lagrange_timeVar)):
    index_timeVar := hasElement(w,lagrange_timeVar):
    if index_timeVar > NumElems(w) then
      w(index_timeVar) := lagrange_timeVar:
      V(i,index_timeVar) := lagrange_constant:
    elif i > Size(V,1) then
      V(i,index_timeVar) := lagrange_constant:
    else
      V(i,index_timeVar) := V(i,index_timeVar) + lagrange_constant:
    end if:
  end do:
end do:
# Substitutierte sinus/cosinus-Terme der Kinematikparameter wieder zurückersetzen

(* Nicht notwendig, wenn oben auskommentiert.
for j from 1 to RowDimension(pkin) do
  V := subs(parse(sprintf("sin_%s", pkin(j,1)))=parse(sprintf("sin(%s)", pkin(j,1))), V):
  V := subs(parse(sprintf("cos_%s", pkin(j,1)))=parse(sprintf("cos(%s)", pkin(j,1))), V):
end do:
*)
#
tmp_tt2 := time():
printf("%s. Komponenten-Matrix des Energie-Regressors berechnet (%dx%d). Rechenzeit: %1.1fs\n", \
  FormatTime("%Y-%m-%d %H:%M:%S"), RowDimension(V), ColumnDimension(V), tmp_tt2-tmp_tt1):
save V, sprintf("../codeexport/%s/tmp/minimal_parameter_calculation_fixb_maple_linsolve_debug1.m", robot_name):
# Rang der Matrix prüfen für Dimension des Parametervektors
tmp_t1 := time():
RankV := Rank(V):
tmp_t2 := time():
printf("%s. Die Komponenten-Matrix des Energie-Regressors hat Rang %d. Rechenzeit: %1.1fs\n", \
  FormatTime("%Y-%m-%d %H:%M:%S"), RankV, tmp_t2-tmp_t1):
# [Diekmeyer2018_S678] Gl. 3.27
tmp_t1 := time():
U := LinearSolve(Transpose(V),Transpose(V)):
tmp_t2 := time():
printf("%s. Basis der Komponenten-Matrix mit LinearSolve gelöst. Resultat ist %dx%d. Rechenzeit: %1.1fs\n", \
  FormatTime("%Y-%m-%d %H:%M:%S"), RowDimension(U), ColumnDimension(U), tmp_t2-tmp_t1):
# Freie Variablen ("_t") entfernen
for i from 1 to Size(V,1) - RankV do
  for j from 1 to Size(V,1) do
    # printf("subs %d,%d\n", i,j):
    U := subs([_t[i, j]=0],U):
  end do:
end do:
# [Diekmeyer2018_S678] Gl. 3.29
Pb, Pd :=permuteEmptyRows(U):
P := Matrix([[Pb],[Pd]]):
save P, Pb, Pd, U, V, RankV, sprintf("../codeexport/%s/tmp/minimal_parameter_calculation_fixb_maple_linsolve_debug2.m", robot_name):
read sprintf("../codeexport/%s/tmp/minimal_parameter_calculation_fixb_maple_linsolve_debug2.m", robot_name):
# [Diekmeyer2018_S678] Gl. 3.31 (links)
Paramvec2 := Pb.U.PV2_vec(11..,..):
# Nachverarbeitung: Prüfe ob Imaginärteile aus Substitutionsalgorithmus vorhanden sind
tmp_t1 := time():
# Diese Prüfung ist nur notwendig, wenn oben die Lagrange-Terme in Exponent-Form gebracht wurden
# Dann ist es möglich, dass imaginäre Ausdrücke im Parametervektor stehen.
imagpart:=false:
for i from 1 to RowDimension(Paramvec2) do
  # Umrechnung der Exponentialform in trigonometrische Form
  Paramvec2(i,1) := simplify(convert(expand(Paramvec2(i,1)), trig));

  if not evalc(Im(Paramvec2[i,1])) = 0 then # müssen eckige Klammern sein, damit evalc funktioniert.
    imagpart := true: # Imaginärteil erkannt. Annahme: Auch nachfolgende Vereinfachung kann daran nichts ändern.
    break:
  end if:
  if has(Paramvec2(i,1), I) then # Prüfe auf Imaginärteil durch Existenz der imaginären Einheit
    Paramvec2(i,1) := evalc(Re(Paramvec2[i,1])): # da Im() oben Null war, muss Re() den vollen Term enthalten.
  end if:
end do:
if imagpart then
  printf("Der Parametervektor hat trotz versuchter Substitution immer noch Imaginärteile aus dem Algorithmus. Verwerfe Lösung und Abbruch.\n"):
  quit: # Funktioniert in GUI nicht richtig...
  robot_name := "": # ...Daher auch Löschung des Roboternamens.
end if:
tmp_t2 := time():
Paramvec2 := simplify2(Paramvec2): # nochmal vereinfachen
;
printf("%s. Parameter-Vektor nachverarbeitet. Rechenzeit: %1.1fs\n", \
  FormatTime("%Y-%m-%d %H:%M:%S"), tmp_t2-tmp_t1):
# Debug: Fertigen Parameter von anderem Durchlauf laden
# read sprintf("../codeexport/%s/tmp/minimal_parameter_vector_fixb_maple_linsolve", robot_name):
printf("%s. Dimension des Minimalparametervektors: %dx%d\n", \
  FormatTime("%Y-%m-%d %H:%M:%S"), RowDimension(Paramvec2), ColumnDimension(Paramvec2)):
Paramvec2;
# [Diekmeyer2018_S678] Gl. 3.31 (rechts)
t_ges_minpar:=Transpose(Pb.Transpose(t_ges)):
u_ges_minpar:=Transpose(Pb.Transpose(u_ges)):
printf("%s. Dimension der Regressormatrix: %dx%d\n", \
  FormatTime("%Y-%m-%d %H:%M:%S"), RowDimension(t_ges_minpar), ColumnDimension(t_ges_minpar)):
# Symbolischer Test des Erfolgs der Parameterminimierung: Ist konstant bzw. Null.
# test_u := simplify2(u_ges_minpar.Paramvec2 - u_ges.PV2_vec(11..,..)):
# test_t := simplify2(t_ges_minpar.Paramvec2 - t_ges.PV2_vec(11..,..)):
# Code exportieren
save Paramvec2, sprintf("../codeexport/%s/tmp/minimal_parameter_vector_fixb_maple", robot_name):
save Paramvec2, sprintf("../codeexport/%s/tmp/minimal_parameter_vector_fixb_maple_linsolve", robot_name): # zum Testen gegen andere Implementierung
if codegen_act then
   MatlabExport(Paramvec2, sprintf("../codeexport/%s/tmp/minimal_parameter_vector_fixb_matlab.m", robot_name), codegen_opt):
end if;
save t_ges_minpar, sprintf("../codeexport/%s/tmp/energy_kinetic_fixb_regressor_minpar_maple.m", robot_name):
save u_ges_minpar, sprintf("../codeexport/%s/tmp/energy_potential_fixb_regressor_minpar_maple.m", robot_name):
if codegen_act then
  MatlabExport(convert_t_s(t_ges_minpar), sprintf("../codeexport/%s/tmp/energy_kinetic_fixb_regressor_minpar_matlab.m", robot_name), codegen_opt):
end if:
if codegen_act then
  MatlabExport(convert_t_s(u_ges_minpar), sprintf("../codeexport/%s/tmp/energy_potential_fixb_regressor_minpar_matlab.m", robot_name), codegen_opt):
end if:

max_stars_repo_name: SchapplM/robsynth-modelgen

max_issues_repo_name: SchapplM/robsynth-modelgen

max_forks_repo_name: SchapplM/robsynth-modelgen


text:
$ifndef _CLASSIFYER_
$define _CLASSIFYER_

$include "Basic.mpl"
$include "Condition.mpl"
$include "Fetch.mpl"
$include "InvOrder.mpl"
$include "InvSimplify.mpl"
$include "InvSol.mpl"
$include "Logout.mpl"
$include "Utils.mpl"

ClassifyHolder:=module()
    option object;
    export
        cid,        # 当前不变量方程中的常数项下标
        ieqCode,    # 不变量方程编号
        sols,       # 当前所有解
        usols,      # 上一个全零不变量方程的解
        oldSols;    # 上一次getSols的解
end module:

# 状态重置
reset:=proc()
    ClassifyHolder:-cid:=0;
    ClassifyHolder:-ieqCode:=0;
    ClassifyHolder:-sols:={};
    ClassifyHolder:-usols:=table();
    ClassifyHolder:-oldSols:={};
end proc:

# 新增解
addSol:=proc(sol)
    ClassifyHolder:-sols:=ClassifyHolder:-sols union {sol};
    return;
end proc:

# 获取新增的代表元
# 新的代表元只能获取一次
getNewSols:=proc()
    local res;
    res:=sort([ (ClassifyHolder:-sols minus ClassifyHolder:-oldSols)[] ],'key'=(x->x:-ieqCode));
    ClassifyHolder:-oldSols:=ClassifyHolder:-sols;
    return res;
end proc:

getCname:=proc()
    ClassifyHolder:-cid:=ClassifyHolder:-cid+1;
    return c[ClassifyHolder:-cid];
end proc:

getIeqCode:=proc()
    ClassifyHolder:-ieqCode:=ClassifyHolder:-ieqCode+1;
    return ClassifyHolder:-ieqCode;
end proc:

classify:=proc(A,As,eqs)
    local sol;
    reset();
    sol:=Object(InvSol):
    sol:-stateCode:=1:
    sol:-oeq:=eqs:
    sol:-As:=As:
    sol:-A:=A:
    sol:-nvars:=LinearAlgebra[RowDimension](A):
    sol:-vars:=[seq(a[i],i=1..sol:-nvars)]:
    resolve(sol);
    return;
end proc:

# 暂时没做重复代表元的处理
getSols:=proc()
    ClassifyHolder:-oldSols:=ClassifyHolder:-sols;
    return sort([ClassifyHolder:-sols[]],'key'=(x->[x:-ieqCode,convert(getDesc(x),`global`)]));
end proc:

resolve:=proc(sol::InvSol)
    local spos,pos,nDelta,_usols,_usol,oldDeltas,oldSol;

    if (sol:-stateCode=1) then
        # 尝试求解偏微分方程组
        # 如果所有方程组为空，则停止求解
        if (sol:-oeq={}) then
            return;
        end if;
        nDelta:=getInvariants(sol:-oeq);
        # 求解失败
        if (indets(nDelta,name) intersect {seq(a[i],i=1..sol:-nvars)} = {}) then
            # 求解失败不添加解
            # 不考虑不能求解不变量的情况
            return "新不变量求解失败";
        end if;
        # 设置新的不变量
        spos:=numelems(sol:-Delta)+1;
        oldDeltas:={sol:-Delta[]};
        oldSol:=Object(sol);
        if (sol:-Delta<>[]) then
            sol:-Delta:=[sol:-Delta[],nDelta[]]:
            # 整体化简不变量
            if (1>=LogLevelHolder:-logLevel) then
                flogf[1]("-----------------------------------------------");
                flogf[1]("对新增不变量按照原不变量进行化简");
                flogf[1]("化简前");
                printDeltas(sol:-Delta);
                sol:-Delta:=simplifyInvariants(sol:-Delta);
                flogf[1]("化简后");
                printDeltas(sol:-Delta);
            end if;
        else
            sol:-Delta:=[sol:-Delta[],nDelta[]]:
        end if;
        sol:-orders:=findInvariantsOrder~(sol:-Delta):# 计算不变量阶数
        # 根据新的不变量是否有约束来决定是否求解上一个全零方程
        nDelta:=remove(x->(x in oldDeltas),sol:-Delta);
        if ( `union`(findDomain~(nDelta)[]) <> {} ) and (sol:-Delta<>[]) then
            solveRestAllZeroIeqs(oldSol);
        end if;
        # 建立和求解不变量方程组
        for pos from spos to numelems(sol:-Delta) do
            buildInvEqs(sol,pos);
        end do;
        # 生成新的不变量
        genInvariants(sol);
    elif (sol:-stateCode=2) then
        # 求解不变量方程组
        solveInvEqs(sol);
    elif (sol:-stateCode=3) then
        # 取代表元
        fetchRep(sol);
    elif (sol:-stateCode=4) then
        # 求解变换方程
        solveTransEq(sol);
    end if;
    return;
end proc:

# 建立不变量的方程组
buildInvEqs:=proc(_sol::InvSol,pos::posint)
    global sols,cid;
    local sol,rs,i,n,x,xpos,eqs;
    n:=numelems(_sol:-Delta);
    # 分奇偶讨论
    if type((_sol:-orders[pos]),even) then
        xpos:=[1,-1,0];
    else
        xpos:=[1,0];
    end if;
    # 生成方程右端
    cid:=0;
    rs:=Array(1..n,x->
    if (x>pos) then
        getCname()
    else
        0
    end if);
    # 逐个方程求解
    for x in xpos do
        # 对于Delta[pos]=0，构建下一个方程进行求解
        # 不求解全零方程
        if (x=0) then
            # 这里是每个全零方程都进行求解的意思
            # 否则直接next就好了
            if (pos<>n) then
                next;
            else
                # 这里是直接求解
                # solveAllZero(_sol);

                # 这里是延后求解
                ClassifyHolder:-usols[getUsolsKey(_sol)]:=_sol;
                return;
            end if;

            # next;
        end if;
        rs[pos]:=x;
        eqs:=[seq(_sol:-Delta[i]=rs[i],i=1..n)];
        sol:=Object(_sol);
        sol:-ieqCode:=getIeqCode();
        sol:-ieq:=eqs;
        sol:-stateCode:=2;
        resolve(sol);
    end do;
    return;
end proc:

# 生成新的代表元
genInvariants:=proc(_sol::InvSol)
    local isols,isol,sol;
    sol:=Object(_sol);
    isols:=ieqsolve(sol:-Delta,[seq(a[i],i=1..sol:-nvars)]);
    # 全部求解失败，则求解上一个全零方程
    if andmap(isol->(subsOeq(sol,isol)="新不变量求解失败"),isols) then
        solveRestAllZeroIeqs(sol);
    end if;
end proc:

# 生成新的不变量方程
# 这么写会导致和非自由变量有关的偏导都变成0
subsOeq:=proc(_sol::InvSol,isol)
    local oeq,sol,v,vv,vars,Delta;
    flogf[1]("--------------------------------------------------------------");
    flogf[1]("求解新的不变量");
    flog[1]({seq(Delta[i]=0,i=1..numelems(_sol:-Delta))});
    flog[2](getDisplayDelta(_sol));
    flogf[1]("取解");
    flog[1](isol);
    oeq:=_sol:-oeq;
    vars:=_sol:-vars;
    v,vv:=selectremove(x->(lhs(x)<>rhs(x)),isol);
    vv:=lhs~(vv);# 方程中的剩余自由变量
    oeq:=PDETools:-dsubs(phi(vars[])=phi(vv[]),oeq);
    oeq:=eval(subs(v[],oeq)) minus {0};
    sol:=Object(_sol);
    sol:-oisol:=isol;
    sol:-stateCode:=1;
    sol:-oeq:=oeq;
    sol:-vars:=vv;
    return resolve(sol);
end proc:

# 求解不变量方程组
solveInvEqs:=proc(_sol::InvSol)
    local isols,icons,n,vars,sol,i;
    n:=_sol:-nvars;
    vars:=[seq(a[i],i=1..n)];
    isols:=ieqsolve(_sol:-ieq,vars);
    icons:=findSolutionDomain~(isols);
    n:=numelems(isols);
    for i from 1 to n do
        sol:=Object(_sol);
        sol:-stateCode:=3;
        sol:-isol:=isols[i];
        sol:-icon:=icons[i];
        resolve(sol);
    end do;
    return;
end proc:

# 对不变量全为0的方程进行求解
solveAllZero:=proc(_sol)
    local sol,var,isols,icons,i,n,reps,rep,nsol,nnsol;
    sol:=Object(_sol);
    sol:-ieq:=[seq(x=0,x in sol:-Delta)];
    sol:-ieqCode:=getIeqCode();
    var:=[seq(a[i],i=1..sol:-nvars)];
    isols:=ieqsolve(sol:-Delta,var);
    icons:=findSolutionDomain~(isols);
    n:=numelems(isols);
    for i from 1 to n do
        nsol:=Object(sol);
        nsol:-isol:=isols[i];
        nsol:-icon:=icons[i];
        reps:=fetchSolRep(nsol,nonzero);
        for rep in reps do
            nnsol:=Object(nsol);
            nnsol:-stateCode:=4;
            setRep(nnsol,rep);
            flogf[1]("--------------------------------------------------------------");
            flogf[1]("求解全零方程");
            flog[1](getDisplayIeq(nnsol));
            flog[1](getDisplayDelta(nnsol));
            flogf[1]("取解");
            flog[1](nnsol:-isol);
            flogf[1]("具有约束条件");
            flog[1](nnsol:-icon);
            flogf[1]("取特解");
            flog[1](nnsol:-rvec);
            flogf[1]("取代表元");
            flog[1](nnsol:-rep);
            resolve(nnsol);
        end do;
    end do;
end proc:

# 取代表元
fetchRep:=proc(_sol::InvSol)
    local n,_ax;
    flogf[1]("--------------------------------------------------------------");
    flogf[1]("对于不变量方程");
    flog[1](getDisplayIeq(_sol));
    flog[1](getDisplayDelta(_sol));
    flogf[1]("取解");
    flog[1](_sol:-isol);
    flogf[1]("具有约束条件");
    flog[1](_sol:-icon);
    n:=_sol:-nvars;
    _ax:=fetchSolRep(_sol);
    if (_ax=NULL) then# 取特解失败
        addSol(_sol);
        flogf[1]("取特解失败");
        return;
    end if;
    setRep(_sol,_ax);
    if (_sol:-rep=0) then
        flogf[1]("代表元取0");
        return;
    end if;
    _ax:=Matrix(_ax);
    _sol:-stateCode:=4;
    flogf[1]("取特解");
    flog[1](convert(_ax,list));
    flogf[1]("取代表元");
    flog[1](_sol:-rep);
    resolve(_sol);
end proc:

solveTransEq:=proc(_sol::InvSol)
    local ax,_ax,n,eq,sol,con;
    n:=_sol:-nvars;
    ax:=Matrix([seq(a[i],i=1..n)]);
    _ax:=_sol:-rvec;
    # a_=a.A
    _sol:-teq[1],_sol:-tsol[1],_sol:-tcon[1]:=solveTeq(_ax,ax,_sol);
    # a=a_.A
    _sol:-teq[2],_sol:-tsol[2],_sol:-tcon[2]:=solveTeq(ax,_ax,_sol);
    if andmap(x->(x=[]),_sol:-tsol) then
        # 无解
        flogf[1]("变换方程求解失败");
        addSol(_sol);
    else
        # 有解
        flogf[1]("变换方程有解");
        _sol:-stateCode:=5;
        addSol(_sol);
        # 在logLevel为1时输出
        if (1>=LogLevelHolder:-logLevel) then
            printTeq(_sol,1);
            printTeq(_sol,2);
        end if;
    end if;
    return;
end proc:

(*
    求解变换方程

    ？ 因为变换方程只关心其存在性，而不在乎其完整性，因此不用explicit选项，而采用convert/radical

    是否选择explicit选项也是有待衡量的，example3出现了不用explicit选项则有复杂约束不能简单消去的情况。
    但是选择explicit选项之后，会面临解过多的问题
    两害取其轻，还是取吧

    该二次求解算法的求解能力值得肯定，测试是否可以转化也依赖于该算法，最好不要修改。
    当然也出现过解的的解不满足方程的情况，至于是否要加以验证，还有待考虑。
*)
solveTeq:=proc(a,b,sol)
    local var,teq,tsol,tcon,scon,eqs,eq,_eq,_con,_sol;
    teq:=convert((a-b.sol:-A),list);
    teq:=subs(sol:-isol[],teq);
    var:=[seq(epsilon[i],i=1..sol:-nvars)];
    tsol:=teqsolve(teq,var);
    if (tsol=[]) then
        # 求解失败，尝试二次求解法方法
        # 首次求解
        eqs:=teqsolve(teq);
        # 二次求解
        tsol:=[];
        tcon:=[];
        for eq in eqs do
            _eq:=select(eqOfEpsilon,eq);
            _con:=remove(eqOfEpsilon,eq);
            _con:=remove(x->type(x,`=`) and (lhs(x)=rhs(x)),_con);
            _sol:=teqsolve(_eq,var,_explicit);
            _con:=map(x->clearConditions(findSolutionDomain(x)) union _con,_sol);
            tsol:=[tsol[],_sol[]];
            tcon:=[tcon[],_con[]];
        end do;
    else
        # 求解成功，直接计算约束
        tcon:=map(x->clearConditions(findSolutionDomain(x)),tsol);
    end if;
    # 清理矛盾解
    tsol:=zip((s,c)->if (undefined in rhs~(c)) then NULL else s end if,tsol,tcon);
    tcon:=remove(c->(undefined in rhs~(c)),tcon);
    return teq,tsol,tcon;
end proc:

eqOfEpsilon:=proc(eq)
    return ormap(x->type(x,specindex(epsilon)),indets(eq,name));
end proc:

# 保留和a,c有关的约束
clearConditions:=proc(con)
    return select(has,con,{a,c});
end proc:

# 求解剩余全零不变量方程
solveRestAllZeroIeqs:=proc(sol::InvSol)
    solveAllZero( ClassifyHolder:-usols[getUsolsKey(sol)] );
end proc:

# 获取全零不变量方程在usols中的key
getUsolsKey:=proc(sol)
    return convert([seq(Delta[i],i=1..numelems(sol:-Delta))],`global`);
end proc:

# 自定义不变量方程求解函数
#
# 采用set的方式指定变量相比于list方式指定变量，解可能更简洁。
# 以为list方式指定变量会优先用后面的变量表示前面的变量，有时这是不好的。
#
# 采用set求解时，example3出现了没有导出a[4]=0的解的情况，例子并不完整。
# 采用list求解时，example包含了a[3]a[4]的约束，也没求解完整。
#
# 可能可以通过自定义补全解的算法来解决这个问题。
# 不过目前只能想到处理不等于0的约束，至于大于0小于0的约束则还没想法。
# 不等于0的约束，可以直接取为0之后再进行求解
ieqsolve:=proc(eq,vars)
    # return convert~([RealDomain:-solve(eq,convert(vars,set),explicit)],list);
    return RealDomain:-solve(eq,vars,explicit);
end proc:

# 自定义变换方程求解函数
teqsolve:=proc({_explicit::boolean:=false})
    if _explicit then
        if _nrest=1 then
            return [RealDomain:-solve(_rest,explicit)];
        elif _nrest=2 then
            return RealDomain:-solve(_rest,explicit);
        else
            error "未知调用方式";
        end if;
    else
        if _nrest=1 then
            return convert~([RealDomain:-solve(_rest)],radical);
        elif _nrest=2 then
            return convert~(RealDomain:-solve(_rest),radical);
        else
            error "未知调用方式";
        end if;
    end if;
end proc:

# 输出Delta
printDeltas:=proc(ds)
    map(i->print(Delta[i]=ds[i]),[seq(x,x=1..numelems(ds))]);
end proc:

$endif


text:
(*
 Copyright (C) 2017 M.A.L. Marques

 This Source Code Form is subject to the terms of the Mozilla Public
 License, v. 2.0. If a copy of the MPL was not distributed with this
 file, You can obtain one at http://mozilla.org/MPL/2.0/.
*)

(* type: gga_exc *)

$define xc_dimensions_2d

_2d_b88_beta := 0.007:
_2d_b88_csi  := 8: (* for harmonic potentials *)

_2d_b88_f := x -> 1 + _2d_b88_beta/X_FACTOR_2D_C*x^2/(1 + _2d_b88_csi*_2d_b88_beta*x*arcsinh(x)):

f := (rs, zeta, xt, xs0, xs1) -> gga_exchange(_2d_b88_f, rs, zeta, xs0, xs1):


text:
`is_element/tree_FFbar` := (N::posint) -> (A::set) -> (TT) -> proc(Q)
 local TT1,B,C,P,T,U;
 global reason;

 if not(type(Q,table)) then
  reason := [convert(procname,string),"Q is not a table",Q];
  return false;
 fi;

 P := sort(map(sort,map(op,{indices(Q)})));
 TT1 := select(T -> nops(T) > 1,TT);

 if P <> TT1 then
  reason := [convert(procname,string),"Q is not indexed by the big sets in TT",P,TT1];
  return false;
 fi;

 C := children_map(A)(TT);

 for T in P do
  if not(`is_element/SCP`(N)(T)(Q[T])) then
   reason := [convert(procname,string),"Q[T] is not in SCP(N)(T)",eval(Q[T]),N,T,reason];
   return false;
  fi;

  for U in C[T] do
   if nops(U) > 1 and not(`is_element/SCP2`(N)(T,U)([Q[T],Q[U]])) then
    reason := [convert(procname,string),"(Q[T],Q[U]) is not in SCP2(N)(T,U)",eval(Q[T]),eval(Q[U]),N,T,U,reason];
    return false;
   fi;
  od;
 od;

 return true;
end;

######################################################################

`random_element/tree_FFbar` := (N::posint) -> (A::set) -> (TT) -> proc()
 `mu/tree_Fbar/tree_FFbar`(N)(A)(TT)(`random_element/tree_Fbar`(N)(A)(TT)());
end;

######################################################################

`is_equal/tree_FFbar` := (N::posint) -> (A::set) -> (TT) -> proc(P,Q)
 local T,TT1;
 global reason;

 TT1 := select(T -> nops(T) > 1,TT);

 for T in TT1 do
  if not(`is_equal/SCP`(N)(T)(P[T],Q[T])) then
   reason := [convert(procname,string),"P[T] <> Q[T]",T,P[T],Q[T],reason];
   return false;
  fi;
 od;

 return true;
end;

######################################################################

`is_leq/tree_FFbar` := (N::posint) -> (A::set) -> (TT) -> proc(P,Q)
 local T,TT1;
 global reason;

 TT1 := select(T -> nops(T) > 1,TT);

 for T in TT1 do
  if not(`is_leq/SCP`(N)(T)(P[T],Q[T])) then
   reason := [convert(procname,string),"P[T] is not <= Q[T]",T,P[T],Q[T],reason];
   return false;
  fi;
 od;

 return true;
end;

######################################################################

`is_interior/tree_FFbar` := (N::posint) -> (A::set) -> (TT) -> proc(Q)
 return `is_element/ICP`(N)(A)(Q[A]);
end;

######################################################################

`inc/ICP/tree_FFbar` := (N::posint) -> (A::set) -> (TT) -> proc(Q0)
 local Q,TT1,T,t;

 Q := table();
 TT1 := select(T -> nops(T) > 1,TT);
 for T in TT1 do
  Q[T] := `top/autorel`(T) intersect Q0;
 od;

 return eval(Q);
end:

######################################################################

`res/tree_FFbar/ICP` := (N::posint) -> (A::set) -> (TT) -> proc(Q)
 if not(`is_interior/tree_FFbar`(N)(A)(TT)(Q)) then
  return FAIL;
 fi;

 return Q[A];
end;

######################################################################

`res/FFbar/tree_FFbar` := (N::posint) -> (A::set) -> (TT) -> proc(Q)
 local TT1,Q1,T;

 TT1 := select(T -> nops(T) > 1,TT);

 Q1 := table():
 for T in TT1 do Q1[T] := eval(Q[T]); od;

 return eval(Q1);
end:

######################################################################

`ext/tree_FFbar/FFbar` := (N::posint) -> (A::set) -> (TT) -> proc(Q)
 local TT1,UU,TU,T,U,m,Q1,i;

 TT1 := select(T -> nops(T) > 1,TT);
 UU := `list_elements/big_subsets`(A);

 Q1 := table():

 for U in UU do
  TU := select(T -> U minus T = {},TT1);
  m := min(op(map(nops,TU)));
  TU := select(U -> nops(U) = m,TU);
  T := TU[1];
  Q1[U] := [seq(`top/autorel`(U) intersect Q[T][i],i=1..N)];
  if not(`is_element/SCP`(N)(U)(Q1[U])) then
   return FAIL;
  fi;
 od;

 return eval(Q1);
end:

######################################################################

`C/tree_FFbar` := (N::posint) -> (A::set) -> (TT) -> proc(Q)
 local L,TT1,T,G,i,ab;
 L := NULL;
 TT1 := select(T -> nops(T) > 1,TT);
 for T in TT1 do
  G := `gamma/SCP`(N)(T)(Q[T]);
  L := L,seq(seq([i,op(ab)],ab in G[i]),i=1..N);
 od;
 return {L};
end:

######################################################################

`list_elements/tree_FFbar` := NULL;
`count_elements/tree_FFbar` := NULL;

######################################################################

`mu/tree_Fbar/tree_FFbar` := (N::posint) -> (A::set) -> (TT) -> proc(x)
 local Q,TT1,T;

 TT1 := select(T -> nops(T) > 1,TT);

 Q := table();

 for T in TT1 do
  Q[T] := `mu/W/ACP`(N)(T)(x[T]);
 od;

 return eval(Q);
end;

######################################################################

`sigma/tree_FFbar/tree_Fbar` := (N::posint) -> (A::set) -> (TT) -> proc(Q)
 local x,TT1,T;

 TT1 := select(T -> nops(T) > 1,TT);

 x := table();

 for T in TT1 do
  x[T] := `sigma/ACP/W`(N)(T)(Q[T]);
 od;

 return eval(x);
end;


text:
(*
 Copyright (C) 2017 M.A.L. Marques

 This Source Code Form is subject to the terms of the Mozilla Public
 License, v. 2.0. If a copy of the MPL was not distributed with this
 file, You can obtain one at http://mozilla.org/MPL/2.0/.
*)

tpss_csi2 := (z, xt, xs0, xs1) ->
  (1 - z^2)*(t_total(z, xs0^2, xs1^2) - xt^2)/(2*(3*Pi^2)^(1/3))^2:

tpss_C0 := (cc, z, xt, xs0, xs1) ->
  + add(cc[i]*z^(2*(i-1)), i=1..4)
  / (1 + tpss_csi2(z, xt, xs0, xs1)*((1 + z)^(-4/3) + (1 - z)^(-4/3))/2)^4:

tpss_aux := (z, xt, ts0, ts1) ->
  xt^2/(8*t_total(z, ts0, ts1)):

tpss_perp := (f_gga, rs, z, xt, xs0, xs1, ts0, ts1) ->
  (1 + tpss_C0(params_a_C0_c, z, xt, xs0, xs1)*tpss_aux(z, xt, ts0, ts1)^2)
  * f_gga(rs, z, xt, xs0, xs1):

tpss_par  := (f_gga, rs, z, xt, xs0, xs1, ts0, ts1) ->
  - (1 + tpss_C0(params_a_C0_c, z, xt, xs0, xs1))*tpss_aux(z, xt, ts0, ts1)^2*(
    + m_max(gga_stoll_par(f_gga, rs,  z, xs0,  1), f_gga(rs, z, xt, xs0, xs1))
    + m_max(gga_stoll_par(f_gga, rs, -z, xs1, -1), f_gga(rs, z, xt, xs0, xs1))
  ):

f_tpss0 := (f_gga, rs, z, xt, xs0, xs1, ts0, ts1) ->
  + tpss_par (f_gga, rs, z, xt, xs0, xs1, ts0, ts1)
  + tpss_perp(f_gga, rs, z, xt, xs0, xs1, ts0, ts1):

f_tpss := (f_gga, rs, z, xt, xs0, xs1, ts0, ts1) ->
  + f_tpss0(f_gga, rs, z, xt, xs0, xs1, ts0, ts1)
  * (1 + params_a_d*f_tpss0(f_gga, rs, z, xt, xs0, xs1, ts0, ts1)*tpss_aux(z, xt, ts0, ts1)^3)
:

max_stars_repo_name: rdietric/lsms

max_issues_repo_name: rdietric/lsms

max_forks_repo_name: rdietric/lsms


text:
(*
 Copyright (C) 2017 M.A.L. Marques

 This Source Code Form is subject to the terms of the Mozilla Public
 License, v. 2.0. If a copy of the MPL was not distributed with this
 file, You can obtain one at http://mozilla.org/MPL/2.0/.
*)

(* type: work_gga_x *)

a1 :=   0.041106:
a2 :=   2.626712:
a3 :=   0.092070:
a4 :=   0.657946:

f0 := s-> a1 * s^a2/(1 + a3 * s^a2)^a4:
f  := x-> f0(X2S*x):


text:
# Definitions for Robot Dynamics Code Generation
# Init
# Erstelle Definitionen für die Maple-Skripte zur Berechnung von Kinematik und Dynamik des Roboters
# Nutze Euler-Winkel (RPY-Konvention) zur Darstellung der Basisorientierung
#
# Quellen:
# [1] Sousa, C. D. and Cortesao, R.: Physical feasibility of robot base inertial parameter identification: A linear matrix inequality approach (2014)
# [2] Ayusawa, K. and Venture, G. and Nakamura, Y.: Identifiability and identification of inertial parameters using the underactuated base-link dynamics for legged multibody systems (2013)
# [3] Fujimoto, Y. and Obata, S. and Kawamura, A.: Robust biped walking with active interaction control between foot and ground (1998)
# [4] Khalil, W. and Kleinfinger, J.-F.: Minimum operations and minimum parameters of the dynamic models of tree structure robots (1987)
# [5] Ortmaier: Robotik I Vorlesungsskript
#
# Moritz Schappler, schappler@irt.uni-hannover.de, 2016-04
# (C) Institut fuer Regelungstechnik, Leibniz Universitaet Hannover
interface(warnlevel=0): # Unterdrücke die folgende Warnung.
restart: # Gibt eine Warnung, wenn über Terminal-Maple mit read gestartet wird.
interface(warnlevel=3):
# Funktion für Euler-Transformation aus Robotik-Repo laden (Pfad muss bekannt sein)
read("../robotics_repo_path"):
read(sprintf("%s/transformation/maple/proc_eulxyzjac", robotics_repo_path)):
# Lese Umgebungsvariable für Codegenerierung.
read "../robot_codegen_definitions/robot_env":
printf("Generiere Parameter für %s\n",robot_name):
# Lese allgemeine Definitionen für Floating-Base
read sprintf("../codeexport/%s/tmp/tree_floatb_twist_definitions", robot_name):
# Robotics Definitions
# Position und Orientierung der Basis (pelvis). Die Orientierung ist mit XYZ-Euler-Winkeln definiert.
# Eine Invertierung der Orientierungsdarstellung sollte nicht notwendig werden, von daher kein Problem mit Orientierungsrepräsentationssingularität.
# gem. [2], S. 4 X_base_t SE(3): Position und Orientierung.
NQB := 6:
X_base_t:=Matrix(6, 1, [rx_base(t),ry_base(t),rz_base(t),alphax_base(t), betay_base(t), gammaz_base(t)]):
X_base_s:=Matrix(6, 1, [rxs_base,rys_base,rzs_base,alphaxs_base, betays_base, gammazs_base]):
# Geschwindigkeit der Basis
# Twist: Verallgemeinerte Geschwindigkeit
V_base_t:=diff~(X_base_t, t):
V_base_s:=Matrix(6, 1, [vxs_base,vys_base,vzs_base,alphaDx_base, betaDy_base, gammaDz_base]):
# Zeitableitung davon
VD_base_t:=diff~(V_base_t, t):
VD_base_s:=Matrix(6, 1, [vDxs_base,vDys_base,vDzs_base,alphaDDx_base, betaDDy_base, gammaDDz_base]):
# Name der Methode für die Orientierungsdarstellung der Basis
base_method_name := "eulxyz":
# Umrechnung von der Ableitung der Basis-Orientierung zu Winkelgeschwindigkeiten:
# Nutze die Euler-XYZ-Konvention ("RPY").
# Siehe [5], Gl. (4.23)
T_basevel := eulxyzjac(X_base_t(4..6,1)):
# Verallgemeinerte Koordinaten, gem [2], S. 4, [3], S.1
NQ:=NQJ+NQB:
q_t := Matrix(NQ,1, <X_base_t, qJ_t>):
q_s := Matrix(NQ,1, <X_base_s, qJ_s>):
qD_t:= Matrix(NQ,1, <V_base_t, qJD_t>):
qD_s:= Matrix(NQ,1, <V_base_s, qJD_s>):
qDD_t:= Matrix(NQ,1, <VD_base_t, qJDD_t>):
qDD_s:= Matrix(NQ,1, <VD_base_s, qJDD_s>):
# MDH-Gelenkwinkel neu speichern (Definition der verallg. Koordinaten war dort noch nicht bekannt
theta := value(theta):
# Export
save q_t, q_s, qD_t, qD_s, qDD_t, qDD_s, qJ_t, qJ_s, qJD_t, qJD_s, qJDD_t, qJDD_s, g_world, X_base_t, X_base_s, V_base_t, V_base_s, VD_base_t, VD_base_s, qoffset, theta, alpha, d, a,v,b,beta, sigma,mu, M, r_i_i_Si, mr_i_i_Si, I_i_i, I_i_Si, PV2_vec, PV2_mat, robot_name, NQ,NQB,NQJ,NJ,NL, base_method_name, T_basevel, kintmp_t, kintmp_s, sprintf("../codeexport/%s/tmp/tree_floatb_eulxyz_definitions", robot_name):
save q_t, q_s, qD_t, qD_s, qDD_t, qDD_s, qJ_t, qJ_s, qJD_t, qJD_s, qJDD_t, qJDD_s, g_world, X_base_t, X_base_s, V_base_t, V_base_s, VD_base_t, VD_base_s, qoffset, theta, alpha, d, a,v,b,beta, sigma,mu, M, r_i_i_Si, mr_i_i_Si, I_i_i, I_i_Si, PV2_vec, PV2_mat, robot_name, NQ,NQB,NQJ,NJ,NL, base_method_name, T_basevel, kintmp_t, kintmp_s, sprintf("../codeexport/%s/tmp/tree_floatb_definitions", robot_name):


text:
# Generated using bssn_4d_spher_matter.mw
eq_evol_rPIb2:=
diff(rPIb2(t,x) ,t) =
expand( simplify (
(-mass*rpsi(t,x) + 2*(diff(rpsi(t, x), x))/(a1(t, x)^2*x)+(diff(diff(rpsi(t, x), x), x))/a1(t, x)^2-(diff(rpsi(t, x), x))*(diff(a1(t, x), x))/a1(t, x)^3+2*rPIb2(t, x)*beta(t, x)/(alpha(t, x)*b1(t, x)^2*a1(t, x)*x)+beta(t, x)*(diff(rPIb2(t, x), x))/(alpha(t, x)*b1(t, x)^2*a1(t, x))+rPIb2(t, x)*(diff(beta(t, x), x))/(alpha(t, x)*b1(t, x)^2*a1(t, x))+2*(diff(rpsi(t, x), x))*(diff(b1(t, x), x))/(b1(t, x)*a1(t, x)^2)+(diff(rpsi(t, x), x))*(diff(alpha(t, x), x))/(alpha(t, x)*a1(t, x)^2) ) * (alpha(t,x)*b1(t,x)^2*a1(t,x))
) );

eq_evol_iPIb2:=
(diff(iPIb2(t, x), t)) = expand( simplify (

( -mass*ipsi(t,x) + 2*(diff(ipsi(t, x), x))/(a1(t, x)^2*x)+(diff(diff(ipsi(t, x), x), x))/a1(t, x)^2-(diff(ipsi(t, x), x))*(diff(a1(t, x), x))/a1(t, x)^3+2*iPIb2(t, x)*beta(t, x)/(alpha(t, x)*b1(t, x)^2*a1(t, x)*x)+iPIb2(t, x)*(diff(beta(t, x), x))/(alpha(t, x)*b1(t, x)^2*a1(t, x))+2*(diff(ipsi(t, x), x))*(diff(b1(t, x), x))/(b1(t, x)*a1(t, x)^2)+(diff(ipsi(t, x), x))*(diff(alpha(t, x), x))/(alpha(t, x)*a1(t, x)^2)+beta(t, x)*(diff(iPIb2(t, x), x))/(alpha(t, x)*b1(t, x)^2*a1(t, x))) * (alpha(t,x)*b1(t,x)^2*a1(t,x))
) );

eq_evol_rpsi := diff(rpsi(t,x),t) = alpha(t,x)/a1(t,x)*rPI(t,x)+beta(t,x)*rPHI(t,x);

eq_evol_ipsi := diff(ipsi(t,x),t) = alpha(t,x)/a1(t,x)*iPI(t,x)+beta(t,x)*iPHI(t,x);


text:
make_GL3_2 := proc()
 global GL3_2;
 local G,L,E,f,n,i,j,k,g,h;

 G := table():

 L := [[]]:
 for i from 1 to 9 do
  L := [seq(seq([op(u),j],j=0..1),u in L)];
 od:

 # Strange order ensures that 1 comes first, followed by other
 # upper unitriangular matrices.
 f := (u) -> Matrix([[u[7],u[5],u[6]],
                     [u[1],u[8],u[4]],
		     [u[3],u[2],u[9]]]);
 L := select(u -> mods(Determinant(f(u)),2) <> 0,L);
 n := nops(L);
 E := map(f,L);

 G["elements"] := E;

 G["index"] := table():
 for i from 1 to n do G["index"][convert(E[i],listlist)] := i; od:

 G["id"] := IdentityMatrix(3);
 G["o"] := (u,v) -> map(mods,u . v,2);
 G["inv"] := (u) -> map(mods,1 / u,2);

 G["o_table"] := table();
 G["inv_table"] := table();

 for i from 1 to n do
  g := G["elements"][i];
  G["inv_table"][i] := G["index"][convert(G["inv"](g),listlist)];

  for j from 1 to n do
   h := G["elements"][j];
   G["o_table"][i,j] := G["index"][convert(G["o"](g,h),listlist)];
  od:
 od:

 G["ee"] := convert(IdentityMatrix(n),listlist);

 G["ring_id"] := G["ee"][1];

 G["antipode"] := proc(u)
  local n,M,w,i,j;
  n := 168;
  M := eval(GL3_2["inv_table"]);
  w := [seq(u[M[i]],i=1..n)];
  return w;
 end:

 G["convolve"] := proc(u,v)
  local n,J,M,w,i,j,ii;
  n := 168;
  J := eval(GL3_2["inv_table"]);
  M := eval(GL3_2["o_table"]);
  w := [0$n];
  for i from 1 to n do
    if u[i] <> 0 then
    ii := J[i];
    w := w +~ [seq(u[i] * v[M[ii,j]],j=1..n)];
   fi;
  od;

  return w;
 end:

 G["borel_elements"] :=
   [seq(seq(seq(Matrix([[1,i,j],[0,1,k],[0,0,1]]),i=0..1),j=0..1),k=0..1)]:
 G["borel_indices"] :=
   map(b -> G["index"][convert(b,listlist)],G["borel_elements"]);
 G["borel_sum"] := [0$168]:
 for i in G["borel_indices"] do
  G["borel_sum"] := G["borel_sum"] +~ G["ee"][i]:
 od:

 G["coxeter_elements"] := map(permutation_matrix,combinat[permute](3)):
 G["coxeter_indices"] := map(w -> G["index"][convert(w,listlist)],G["coxeter_elements"]):
 G["signed_coxeter_sum"] := [0$168]:
 for i in G["coxeter_indices"] do
  G["signed_coxeter_sum"] := G["signed_coxeter_sum"] +~ Determinant(E[i]) *~ G["ee"][i]:
 od:

 G["steinberg_idempotent"] :=
  G["convolve"](G["signed_coxeter_sum"], G["borel_sum"]) / 21;

 GL3_2 := eval(G);
 return(eval(G));
end:


text:
# Domain types which are registered/unregistered with TypeTools.
# Note that the types have to be quoted (additionally to the quote one would
# normally place on certain types) to work properly
local DomBoundVar_type := proc(nm)
    local ixs := [indices(Domain:-ExtBound, nolist)], vs;
    ormap(i->type(nm,Domain:-ExtBound[i]:-VarType), ixs);
end proc;

local DomBoundRange_type := proc(nm)
    local ixs := [indices(Domain:-ExtBound, nolist)], vs;
    ormap(i->type(nm,Domain:-ExtBound[i]:-RangeType), ixs);
end proc;

local DomainTypes := table(
       # Domain bounds
       [(DomBoundVar = 'thismodule:-DomBoundVar_type')
       ,(DomBoundRange = 'thismodule:-DomBoundRange_type')
       ,(DomBoundBinder = ''DInto(DomBoundVar, DomBoundRange, DomBoundKind)'' )
       ,(DomBoundKind   = 'And(name, satisfies(nm->assigned(Domain:-ExtBound[nm])))' )
       ,(DomBound       = ''And(specfunc(DBound)
                               ,{anyfunc(list(DomBoundBinder))
                                ,anyfunc(list(DomBoundBinder),DomCtx)
                                ,anyfunc(list(DomBoundBinder),DomCtx,table) })'' )
       # Domain shape
       ,(DomConstrain = 'specfunc(relation, `DConstrain`)' )
       ,(DomSum       = 'specfunc(DomShape, `DSum`)' )
       ,(DomSplit     = ''DSplit(Partition(DomShape))'' )
       ,(DomInto      = ''Or(DInto(DomBoundVar, DomBoundRange, DomShape)
                            ,DInto(DomBoundVar, DomShape) )'' )
       ,(DomShape     = 'Or( DomConstrain, DomSum, DomSplit, DomInto )' )
       # Domain
       ,(DomCtx = ''t_kb'')
       ,(HDomain = ''DOMAIN(DomBound, DomShape)'' )
       # Maybe domain
       ,(DomNoSol  = ''Not(freeof(`DNoSol`))'' )
       ,(HDomain_mb = ''Or(HDomain, DOMAIN(DomBound, DomNoSol))'' )
       ] );

max_stars_repo_name: vmchale/hakaru

max_issues_repo_name: vmchale/hakaru

max_forks_repo_name: vmchale/hakaru


text:
local Simplify_DConstrain := (can_simp, do_simp) -> module()
    uses Domain;

    export ModuleApply := proc(vs0 :: DomBound, sh :: DomShape, $)
        local vs := vs0, ctx;
        if _Env_HakaruSolve=false then
           return sh end if;
        vs  := Domain:-Bound:-varsOf(vs,"set");
        ctx := Domain:-Bound:-toConstraints(vs0, 'bound_types');
        subsindets(sh, DomConstrain
                  ,x->do_simpl_constraints(vs0, vs, ctx, x));
    end proc;

    local do_simpl_constraints := proc(vs0, vars, ctx_vs, x, $)
        local ctx1, ctx, ss, td, rest, d, in_vs;
        ss, ctx := selectremove(q->depends(q,indets(vars,Name)), x);
        td, rest := selectremove(can_simp(vs0), ss);
        ctx1 := { op(ctx), op(ctx_vs), op(rest) };
        d := map(x->do_simp(vs0,ctx1,x), td);
        DConstrain(op(d), op(ctx), op(rest));
    end proc;
end module;

constraints_about_vars := module()
  export SimplName  := "Make constraints abouts vars";
  export SimplOrder := 6;
  export ModuleApply := Simplify_DConstrain(can_make_about, try_make_about);

  local can_make_about := proc(vs)
    local vars := Domain:-Bound:-varsOf(vs,"set");
    c -> c::relation and (not(lhs(c) in vars) and not(rhs(c) in vars) and has(c,{ln,exp}));
  end proc;
  local try_make_about := proc(dbnd, ctx1, q0, $)
      local vars_q, vars, q_s, q_r, q := q0;
      vars := Domain:-Bound:-varsOf(dbnd,"set");
      vars_q := indets(q, name) intersect vars;
      if nops(vars_q)=0 then return q end if;

      vars_q := op(1,vars_q);
      q := KB:-try_improve_exp(q, vars_q, ctx1);
      q_r := `if`(has(q,{ln,exp}),q0,q);
      q_s := 'solve({q},[op(vars_q)], 'useassumptions'=true) assuming (op(ctx1))';
      q_s := eval(q_s);
      if q_s::list then
        if nops(q_s)=0 then return q_r end if;
        q_s := map(s->remove(c->c in ctx1 or `and`(c::relation,lhs(c)::name,lhs(c)=rhs(c)), s), q_s);
        q_s := remove(x->x=[],q_s);
        if nops(q_s)=0 then return q_r end if;
        userinfo(3, 'constraints_about_vars'
                ,printf("Found a solution to %a (with solve):\n\t%a\n", q0, op(1,q_s)));
        op(op(1,q_s)); # pick the first solution arbitrarily!
      else
        q_r
      end if;
  end proc;
end module;

clamp_extraneous_constraints := module()
  export SimplName  := "clamp_extraneous_constraints";
  export SimplOrder := 6 - 1/10;
  export ModuleApply := Simplify_DConstrain(can, `try`);

  local can := proc(vs)
    local vars := Domain:-Bound:-varsOf(vs,"set");
    c-> c::relation and
        (ormap(s->s(c)::Name and not(c in vars),[lhs,rhs]) and has(c, vars));
  end proc;

  local `try` := proc(dbnd, ctx1, q0, $)
    local bnd, b_ty, b_rel, b_var, b_bnd, b_bnd_vs, vars, extremum, q := q0;
    if has(q, {ln,exp}) then
      vars := Domain:-Bound:-varsOf(dbnd, "set");
      bnd := Domain:-Improve:-classify_relation(q0, x -> type(x, Name) and not (x in vars));
      if bnd=FAIL then return q0 end if;
      b_ty, b_rel, b_var, b_bnd := op(bnd):
      b_bnd_vs := indets(b_bnd, satisfies(x -> x in vars));
      if b_bnd_vs <> {} then
        b_ty := `if`(b_ty=B_LO, 'minimize', 'maximize');
        extremum := b_ty(b_bnd, op(b_bnd_vs)) assuming(op(ctx1));
        if not (ext :: SymbolicInfinity or has(extremum, {maximize,minimize})) then
          q := q, b_rel(b_var,extremum);
        end if;
      end if;
    end if;
    q
  end proc;
end module;

# Pushes constraints down, or pulls them up, when there are such constraints.
local do_ctx_dir := dir -> proc(vs :: DomBound, sh :: DomShape, $)
    local ctx := `if`(nops(vs)=2,op(2,vs),{});
    ctx := remove(x->x::`::`,ctx);
    if ctx <> {} then
        dir(ctx)(sh);
    else
        sh;
    end if;
end proc;

push_ctx_down := module()
  uses Domain;
  export ModuleApply := do_ctx_dir(ctx->sh->subsindets(sh, DomConstrain, x->if nops(x)<>0 then DConstrain(op(x), op(ctx)) else DConstrain() end if));
  export SimplName  := "Push context down";
  export SimplOrder := 1;
end module;

pull_ctx_out := module()
  export ModuleApply :=
    do_ctx_dir(ctx->sh->subsindets(sh, DomConstrain
              ,x->remove(c->c in ctx or (is(c) assuming op(ctx)), x)));
  export SimplName  := "Pull context out";
  export SimplOrder := 20;
end module;

# Turns constraints(relations) into bounds
classify_DConstrains := module()
  uses Domain;

  export SimplName  := "Classify constraints";
  export SimplOrder := 8;

  export ModuleApply := proc(dbnds::DomBound,dsh::DomShape,$)
    subsindets(dsh, DomConstrain, x->classify_Rels([op(x)], dbnds));
  end proc;

  # Split relations into those we can incorporate
  # into bounds and others; this process also includes
  # some intermediate computations about the new bound:
  #  [ DInto(var), orig, upd_bound, bound_kind, numVars ]
  #    var         :=   variable to be modified
  #    orig        :=   the original relation
  #    upd_bound   :=   how to modify (a function)
  #    bound_kind  :=   {LO/HI}
  #    numVars     :=   how many variables are mentioned by the bound
  local classify_Rel := proc(r0,dbnd,$)
    local ret, vars, mk_k, bnd_k, i, var_k, q, cl, lhs_r, rhs_r, r := r0;
    ret := [ DConstrain, r0 ];
    vars := Domain:-Bound:-varsOf(dbnd,"set");

    # extract a representation of the relation, and check that it is a bound
    cl := Domain:-Improve:-classify_relation(r, vars);
    if cl=FAIL or not (op(1,cl) in {B_LO,B_HI}) then return ret; end if;
    bnd_k,mk_k,lhs_r,rhs_r := op(cl);

    i := Domain:-Bound:-varIx_mb(dbnd, lhs_r);
    if i >= 0 then
      var_k := op([1,i,3], dbnd);
      q := Domain:-ExtBound[var_k]:-RecogBound(mk_k,rhs_r);
      if q <> NULL then
        ret := [ DInto(lhs_r), r0, q, bnd_k, countVs(vars minus {lhs_r},r0) ];
      end if;
    end if;
    ret;
  end proc;

  # Classifies a set of relations as potential bounds. The process
  # is as follows:
  #  - Classify each relation individually (see classify_Rel)
  #  - "Categorize" by the classification
  #  - Split relations which will be DConstrains, and those where we found
  #      more than one LO/HI bound or more than 2 bounds total.
  #  - Sort the relations which are to become bounds by the number of variables
  #      mentioned in the bound
  #  - Rebuild the continuations producing each new bound (see make_DInto)
  #  - Rebuild the leftover relations as a DConstrain.
  #  - Squash together ("foldr (.) id fs z") the continuations and leftover relations.
  local classify_Rels := proc(rs0,dbnd,$)
    local cs1, cs2, cs, rs := rs0;
    rs := map(r->classify_Rel(r,dbnd), rs);
    rs := [ListTools[Categorize]( ((a,b) -> op(1,a)=op(1,b) ), rs )];
    cs1, rs := selectremove(x->op([1,1],x)=DConstrain, rs);
    rs, cs2 := selectremove(x->nops(x)=1 or (nops(x)=2 and op([1,4],x)<>op([2,4],x)), rs);
    rs := sort(rs, key=solOrder);
    rs := map(r->make_DInto(r,dbnd),rs);
    cs := DConstrain(map(c->map(q->op(2,q),c)[], [op(cs1),op(cs2)])[]);
    foldr((f,g)->f@g,x->x,op(rs))(cs);
  end proc;

  # Rebuild a classified (by classify_Rel) relation; this
  # takes a list of such classifications, all of whose variable should
  # be identical, and for which the different bound kinds are all distinct.
  # None of this is checked!
  # (since there are 2 bound kinds, this list is at most length 2)
  local make_DInto := proc(r::list([specfunc(DInto),anything,anything,identical(B_LO,B_HI),nonnegint]), dbnd, $)
    local var_nm, var_ty;
    var_nm := op([1,1,1], r);
    var_ty := Domain:-Bound:-get(dbnd, var_nm)[1];
    var_ty := foldr( ((q,x)->op(3,q)(x)), var_ty, op(r) );
    (ctx -> DInto(var_nm, var_ty, ctx))
  end proc;
  local flip_rel := table([`=`=`=`,`<>`=`<>`,`<=`=`>=`,`<`=`>`,`>=`=`<=`,`>`=`<`]);
  local countVs := ((vs,c)-> nops(indets(c, DomBoundVar) intersect {op(vs)} ));
  local solOrder := (x->`+`(map(z->op(5,z),x)[]));
end module;

# todo; this should actually solve for a variable, then substitute
# that variable in. In most cases, it would probably be enough to
# leave that as it is; it would simplify later.
singular_pts := module()
  export SimplName  := "Single_pts";
  export SimplOrder := 14;

  local can_remove := ((xs,vs_ty) ->
     nops(xs)>0 and
     andmap(x ->
     x :: `=` and
     ormap(s->s(x)::vs_ty,[lhs,rhs]) and
     not(is(x)), xs));

  export ModuleApply := proc(bnds_ :: DomBound, sh_ :: DomShape, $)
    local bnds := bnds_, sh := sh_, vs, todo, sh1, vs_ty;
    vs := applyop(bl -> select(b->op(3,b)=`Int`, bl), 1, bnds);
    vs := Domain:-Bound:-varsOf(vs,"set");
    vs_ty := 'satisfies(x->x in vs)';
    todo := indets(sh, specfunc(`DConstrain`));
    todo := select(x->can_remove(x,vs_ty), todo);
    subs([seq(t=DSum(),t=todo)], sh);
  end proc;
end module;


text:
(*
 Copyright (C) 2017 M.A.L. Marques

 This Source Code Form is subject to the terms of the Mozilla Public
 License, v. 2.0. If a copy of the MPL was not distributed with this
 file, You can obtain one at http://mozilla.org/MPL/2.0/.
*)

(* type: work_mgga_c *)

par_n := 5:

gamma_x := 0.004:
par_x := [
    [  1.000, 0, 0],
    [  1.308, 0, 1],
    [  1.901, 0, 2],
    [  0.416, 1, 0],
    [  3.070, 1, 1]
]:

gamma_ss := 0.2:
par_ss := [
    [  1.000, 0, 0],
    [ -1.855, 0, 2],
    [ -5.668, 1, 0],
    [-20.497, 3, 2],
    [-20.364, 4, 2]
]:

gamma_os := 0.006:
par_os := [
    [  1.000, 0, 0],
    [  1.573, 0, 1],
    [ -6.298, 0, 3],
    [  2.535, 1, 0],
    [ -6.427, 3, 2]
]:

$define lda_x_params
$include "lda_x.mpl"
$include "b97mv.mpl"

f :=  (rs, z, xt, xs0, xs1, ts0, ts1, us0, us1) ->
  f_b97mv(f_lda_x, rs, z, xt, xs0, xs1, ts0, ts1):


text:
(*
 * 保存不变量分类状态的对象
 *
 * 状态分为：
 *     1 线性偏微分方程组求解失败
 *     2 不变量方程组求解失败
 *     3 取特解失败
 *     4 变换方程求解失败
 *     5 变换方程求解成功，特解可能是代表元，是否是代表元留给人工判断
*)
$ifndef _INVSOL_
$define _INVSOL_

InvSol:=module()
    option object;
    export
        # 导出变量
        stateCode::{1,2,3,4,5},     # 状态代码
        oieq:={},                   # 导出新偏微分方程的不变量方程（如果存在的话）
        oisol,                      # 导出新偏微分方程的不变量方程的解（如果存在的话）
        oeq::set,                   # 导出不变量的偏微分方程
        Delta:=[],                  # 不变量
        orders:=[],                 # 不变量对应的阶数
        ieqCode,                    # 不变量方程的代码
        ieq::list,                  # 不变量方程
        isol,                       # 不变量方程的解（一个InvSol对象只含一个解）
        icon::set,                  # 不变量方程解的条件
        teq:=[[],[]],               # 变换方程（有两个）
        tsol:=[[],[]],              # 变换方程的解，每个方程都有可能有多个解
        tcon:=[[],[]],              # 变换方程的解对于的条件
        rep,                        # 代表元的表达式形式
        rvec,                       # 代表元的系数矩阵形式
        vars,                       # 偏微分方程中的剩余自由变量，初始时为a[1]..a[nvars]
        nvars,                      # 生成元的个数
        As::static,                 # 每个生成元的伴随矩阵
        A::static,                  # 总的伴随变换矩阵
        # 导出函数
        ## 输出相关
        getDisplayIeq::static,      # 简化显示的不变量方程，将不变量的实际表达式简写为Delta
        getDisplayDelta::static,    # 输出不变量
        getDesc::static,            # 获取对象的特征描述
        ModulePrint::static,        # 打印和显示对象的值
        printSol::static,           # 详细显示对象信息
        printTeq::static,           # 显示变换方程的结果
        ## 设置相关
        setRep::static,             # 重新取代表元
        setIsol::static,            # 重新对不变量方程取解
        getRep::static;             # 获取代表元，转化为global对象便于比较

    # 重新取代表元
    setRep:=proc(s::InvSol,rvec::list,{nocheck::boolean:=false})
        local v,r;
        if not nocheck then
            # 解的检验
            r:=eval(subs(seq(a[i]=rvec[i],i=1..numelems(rvec)),s:-ieq));
            if not andmap(evalb,r) then
                error "不是不变量方程的解，代入结果为%1",r;
            end if;
        end if;
        s:-stateCode:=4;
        s:-rvec:=Matrix(rvec);
        s:-rep:=add(rvec[i]*v[i],i=1..numelems(rvec));
        return;
    end proc:

    # 重新对不变量方程取解
    setIsol:=proc(s::InvSol,isol,{nocheck::boolean:=false})
        local r;
        if not nocheck then
            # 解的检验
            r:=RealDomain:-simplify(subs(isol[],s:-ieq));
            if not andmap(evalb,r) then
                error "不是不变量方程的解，代入结果为%1",r;
            end if;
        end if;
        s:-stateCode:=3;
        s:-isol:=isol;
        s:-icon:=findSolutionDomain(isol);
        return;
    end proc:

    # 获取代表元，转化为global对象便于比较
    getRep:=proc(s::InvSol)
        return convert(s:-rep,`global`);
    end proc:

    # 简化显示的不变量方程，将不变量的实际表达式简写为Delta
    getDisplayIeq:=proc(self::InvSol)
        local Delta;
        return {seq(Delta[i]=rhs(self:-ieq[i]),i=1..numelems(self:-Delta))};
    end proc:

    # 输出不变量
    getDisplayDelta:=proc(self::InvSol)
        local Delta;
        return [seq(Delta[i]=self:-Delta[i],i=1..numelems(self:-Delta))];
    end proc:

    # 获取对象的特征描述
    getDesc:=proc(s)
        if   (s:-stateCode=1) then
            return s:-oeq;
        elif (s:-stateCode=2) then
            return getDisplayIeq(s);
        elif (s:-stateCode=3) then
            return s:-isol;
        elif (s:-stateCode=4) then
            return s:-rep;
        elif (s:-stateCode=5) then
            return s:-rep;
        end if;
    end proc:

    # 打印和显示对象的值
    ModulePrint:=proc(s)
        return getDesc(s);
    end proc:

    # 输出解
    printSol:=proc(s::InvSol)
        printf("---------------------------------------------------------");
        if     (s:-stateCode=1) then
            printf("新的不变量求解失败，状态代码1");
            print(s:-oieq);
            printf("取解");
            print(s:-oisol);
            printf("求解失败的偏微分方程为");
            print(s:-oeq);
        elif    (s:-stateCode=2) then
            printf("不变量方程求解失败，状态代码2");
            print(getDisplayIeq(s));
        elif    (s:-stateCode=3) then
            printf("取代表元失败，状态代码3");
            print(getDisplayIeq(s));
            printf("取解");
            print(s:-isol);
        elif    (s:-stateCode=4) then
            printf("变换方程求解失败，状态代码4");
            print(getDisplayIeq(s));
            printf("取解");
            print(s:-isol);
            printf("具有约束");
            print(s:-icon);
            printf("取代表元");
            print(s:-rep);
            printf("求解失败的两个变换方程为");
            print~(s:-teq);
        elif    (s:-stateCode=5) then
            printf("变换方程求解成功，状态代码5");
            print(getDisplayIeq(s));
            printf("取解");
            print(s:-isol);
            printf("具有约束");
            print(s:-icon);
            printf("取代表元");
            print(s:-rep);
            printf("变换方程有解");
            printTeq(s,1);
            printTeq(s,2);
        end if;
        printf("---------------------------------------------------------");
        return;
    end proc:

    # 输出变换方程和解
    printTeq:=proc(sol,pos)
        if (sol:-tsol[pos]=[]) then
            printf("变换方程 %d 无解",pos);
            printf("方程为");
            print(sol:-teq[pos]);
        else
            printf("变换方程 %d 有解",pos);
            printf("方程为");
            print(sol:-teq[pos]);
            printf("有解");
            print(sol:-tsol[pos]);
            printf("具有条件");
            print(sol:-tcon[pos]);
        end if;
    end proc:
end module:

$endif


text:
var $arr1 <[10] f32>

func $ArrayAdd () f32 {
 var %ff i32
 var %sum f32
 dassign %sum ( constval f32 0.0f )
 foreachelem %ff $arr1 {
   dassign %sum ( add f32 ( dread f32 %sum, dread f32 %ff ) )
 }
 return ( dread f32 %sum )
}
 # EXEC: %irbuild Main.mpl
 # EXEC: %irbuild Main.irb.mpl
 # EXEC: %cmp Main.irb.mpl Main.irb.irb.mpl

max_stars_repo_name: harmonyos-mirror/OpenArkCompiler-test

max_issues_repo_name: harmonyos-mirror/OpenArkCompiler-test

max_forks_repo_name: harmonyos-mirror/OpenArkCompiler-test


text:
kernelopts(opaquemodules=false):
kernelopts(assertlevel=2):

max_stars_repo_name: lairez/binomsum

max_issues_repo_name: lairez/binomsum

max_forks_repo_name: lairez/binomsum