diff --git a/CHANGELOG.md b/CHANGELOG.md
index 899aa29b..6461a13f 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -4,6 +4,10 @@
 
 * Update memory settings and fix inconsistencies
 
+### Enhancements
+
+* Use multi-line output for PARI/GP formulas
+
 ## v26.3.10
 
 ### Bugfixes
diff --git a/src/form/pari.cpp b/src/form/pari.cpp
index 1f21d575..6f7cc7fd 100644
--- a/src/form/pari.cpp
+++ b/src/form/pari.cpp
@@ -142,7 +142,16 @@ std::string PariFormula::toString() const {
     }
     return buf.str();
   } else {
-    return main_formula.toString("; ", true);
+    // Get all function definitions in the formula to determine output format
+    auto functions =
+        FormulaUtil::getDefinitions(main_formula, Expression::Type::FUNCTION);
+    if (functions.size() > 1) {
+      // Multiple functions: use newline separator and no brackets
+      return main_formula.toString("\n", false);
+    } else {
+      // Single function: no separator needed (single entry)
+      return main_formula.toString("", false);
+    }
   }
 }
 
diff --git a/tests/formula/pari-function.txt b/tests/formula/pari-function.txt
index 178390cb..48fb88e5 100644
--- a/tests/formula/pari-function.txt
+++ b/tests/formula/pari-function.txt
@@ -1,13 +1,22 @@
 A000023: a(n) = if(n==3,-2,if(n==2,2,if(n==1,-1,if(n==0,1,local(l1=a(n-1)); (n-1)*(2*a(n-2)+l1)-l1))))
 A000027: a(n) = n
 A000032: a(n) = if(n==2,3,if(n==1,1,if(n==0,2,a(n-1)+a(n-2))))
-A000058: (a(n) = b(n)+1); (b(n) = if(n==0,1,local(l1=b(n-1)); l1*(l1+1)))
+A000058: a(n) = b(n)+1
+  b(n) = if(n==0,1,local(l1=b(n-1)); l1*(l1+1))
 A000247: a(n) = 2^n-n-2
 A001286: a(n) = floor(((n-1)*n!)/2)
-A001715: (a(n) = truncate(b(n)/6)); (b(n) = if(n==0,1,n*b(n-1)))
+A001715: a(n) = truncate(b(n)/6)
+  b(n) = if(n==0,1,n*b(n-1))
 A001803: a(n) = if((binomial(2*n+1,n)*(n+1))==0,0,(binomial(2*n+1,n)*(n+1))/(2^valuation(binomial(2*n+1,n)*(n+1),2)))
-A002944: (a(n) = b(n-1)); (b(n) = if(n==2,2,if(n==1,1,if(n==0,1,local(l2=b(n-1)); local(l1=c(n-1)); if(gcd(l2,n+1)==0,l1,if((l1%gcd(l2,n+1))==0,l1/gcd(l2,n+1),l1)))))); (c(n) = if(n==2,6,if(n==1,2,if(n==0,1,local(l2=b(n-1)); local(l1=c(n-1)); (n+1)*if(gcd(l2,n+1)==0,l1,if((l1%gcd(l2,n+1))==0,l1/gcd(l2,n+1),l1))))))
+A002944: a(n) = b(n-1)
+  b(n) = if(n==2,2,if(n==1,1,if(n==0,1,local(l2=b(n-1)); local(l1=c(n-1)); if(gcd(l2,n+1)==0,l1,if((l1%gcd(l2,n+1))==0,l1/gcd(l2,n+1),l1)))))
+  c(n) = if(n==2,6,if(n==1,2,if(n==0,1,local(l2=b(n-1)); local(l1=c(n-1)); (n+1)*if(gcd(l2,n+1)==0,l1,if((l1%gcd(l2,n+1))==0,l1/gcd(l2,n+1),l1)))))
 A006430: a(n) = truncate((floor(((-(1==n)+n)*(-(1==n)+n+1)*(-(1==n)+n+2)*(-(1==n)+n+3))/(-(1==n)+n+1))*((-(1==n)+n)*((-(1==n)+n)*((-(1==n)+n)*((-(1==n)+n)*((-(1==n)+n)*((-(1==n)+n)*((-(1==n)+n)*((-(1==n)+n)*((-(1==n)+n)*(23*n-23*(1==n)+963)+17544)+147952)+481675)-1052153)-7850914)-2900162)+60869272)+37067400)-179920800))/79833600)
-A077849: (a(n) = c(n+2)); (b(n) = if(n==3,8,if(n==2,3,if(n==1,1,if(n==0,0,2*b(n-1)+b(n-2)+b(n-3)+d(n-1)))))); (c(n) = if(n==3,3,if(n==2,1,if(n==1,0,if(n==0,0,b(n-1)))))); (d(n) = if(n==3,1,if(n==2,1,if(n==1,1,if(n==0,1,d(n-1))))))
-A215020: (a(n) = c(n-1)); (b(n) = if(n==1,2,if(n==0,1,local(l1=b(n-1)); l1+floor((c(n-1)+1)/(if(l1==0,0,valuation(l1,2))+1))))); (c(n) = if(n==1,0,if(n==0,0,local(l1=b(n-1)); (c(n-1)+1)%(if(l1==0,0,valuation(l1,2))+1))))
+A077849: a(n) = c(n+2)
+  b(n) = if(n==3,8,if(n==2,3,if(n==1,1,if(n==0,0,2*b(n-1)+b(n-2)+b(n-3)+d(n-1)))))
+  c(n) = if(n==3,3,if(n==2,1,if(n==1,0,if(n==0,0,b(n-1)))))
+  d(n) = if(n==3,1,if(n==2,1,if(n==1,1,if(n==0,1,d(n-1)))))
+A215020: a(n) = c(n-1)
+  b(n) = if(n==1,2,if(n==0,1,local(l1=b(n-1)); l1+floor((c(n-1)+1)/(if(l1==0,0,valuation(l1,2))+1))))
+  c(n) = if(n==1,0,if(n==0,0,local(l1=b(n-1)); (c(n-1)+1)%(if(l1==0,0,valuation(l1,2))+1)))
 A254006: a(n) = if(if(((-n)%(-2))==0,(-n)/(-2),-n)<=(-1),0,3^if(((-n)%(-2))==0,(-n)/(-2),-n))