diff --git a/internal/stats/latest_stats.csv b/internal/stats/latest_stats.csv
index eb7b4efb74..f5be61615f 100644
--- a/internal/stats/latest_stats.csv
+++ b/internal/stats/latest_stats.csv
@@ -153,14 +153,14 @@ pairing_bls12377,bls24_315,plonk,0,0
 pairing_bls12377,bls24_317,plonk,0,0
 pairing_bls12377,bw6_761,plonk,51280,51280
 pairing_bls12377,bw6_633,plonk,0,0
-pairing_bls12381,bn254,groth16,1429070,2382640
+pairing_bls12381,bn254,groth16,1419904,2366999
 pairing_bls12381,bls12_377,groth16,0,0
 pairing_bls12381,bls12_381,groth16,0,0
 pairing_bls12381,bls24_315,groth16,0,0
 pairing_bls12381,bls24_317,groth16,0,0
 pairing_bls12381,bw6_761,groth16,0,0
 pairing_bls12381,bw6_633,groth16,0,0
-pairing_bls12381,bn254,plonk,5629807,5285448
+pairing_bls12381,bn254,plonk,5593770,5250897
 pairing_bls12381,bls12_377,plonk,0,0
 pairing_bls12381,bls12_381,plonk,0,0
 pairing_bls12381,bls24_315,plonk,0,0
@@ -181,70 +181,70 @@ pairing_bls24315,bls24_315,plonk,0,0
 pairing_bls24315,bls24_317,plonk,0,0
 pairing_bls24315,bw6_761,plonk,0,0
 pairing_bls24315,bw6_633,plonk,141249,141249
-pairing_bn254,bn254,groth16,969638,1614382
+pairing_bn254,bn254,groth16,963003,1603091
 pairing_bn254,bls12_377,groth16,0,0
 pairing_bn254,bls12_381,groth16,0,0
 pairing_bn254,bls24_315,groth16,0,0
 pairing_bn254,bls24_317,groth16,0,0
 pairing_bn254,bw6_761,groth16,0,0
 pairing_bn254,bw6_633,groth16,0,0
-pairing_bn254,bn254,plonk,3798583,3560759
+pairing_bn254,bn254,plonk,3771397,3534755
 pairing_bn254,bls12_377,plonk,0,0
 pairing_bn254,bls12_381,plonk,0,0
 pairing_bn254,bls24_315,plonk,0,0
 pairing_bn254,bls24_317,plonk,0,0
 pairing_bn254,bw6_761,plonk,0,0
 pairing_bn254,bw6_633,plonk,0,0
-pairing_bw6761,bn254,groth16,3014749,4979960
+pairing_bw6761,bn254,groth16,2592181,4256159
 pairing_bw6761,bls12_377,groth16,0,0
 pairing_bw6761,bls12_381,groth16,0,0
 pairing_bw6761,bls24_315,groth16,0,0
 pairing_bw6761,bls24_317,groth16,0,0
 pairing_bw6761,bw6_761,groth16,0,0
 pairing_bw6761,bw6_633,groth16,0,0
-pairing_bw6761,bn254,plonk,11486969,10777222
+pairing_bw6761,bn254,plonk,9920293,9270827
 pairing_bw6761,bls12_377,plonk,0,0
 pairing_bw6761,bls12_381,plonk,0,0
 pairing_bw6761,bls24_315,plonk,0,0
 pairing_bw6761,bls24_317,plonk,0,0
 pairing_bw6761,bw6_761,plonk,0,0
 pairing_bw6761,bw6_633,plonk,0,0
-scalar_mul_G1_bn254,bn254,groth16,74345,117078
+scalar_mul_G1_bn254,bn254,groth16,69013,108022
 scalar_mul_G1_bn254,bls12_377,groth16,0,0
 scalar_mul_G1_bn254,bls12_381,groth16,0,0
 scalar_mul_G1_bn254,bls24_315,groth16,0,0
 scalar_mul_G1_bn254,bls24_317,groth16,0,0
 scalar_mul_G1_bn254,bw6_761,groth16,0,0
 scalar_mul_G1_bn254,bw6_633,groth16,0,0
-scalar_mul_G1_bn254,bn254,plonk,278909,261995
+scalar_mul_G1_bn254,bn254,plonk,260289,244439
 scalar_mul_G1_bn254,bls12_377,plonk,0,0
 scalar_mul_G1_bn254,bls12_381,plonk,0,0
 scalar_mul_G1_bn254,bls24_315,plonk,0,0
 scalar_mul_G1_bn254,bls24_317,plonk,0,0
 scalar_mul_G1_bn254,bw6_761,plonk,0,0
 scalar_mul_G1_bn254,bw6_633,plonk,0,0
-scalar_mul_P256,bn254,groth16,100828,161106
+scalar_mul_P256,bn254,groth16,93170,148354
 scalar_mul_P256,bls12_377,groth16,0,0
 scalar_mul_P256,bls12_381,groth16,0,0
 scalar_mul_P256,bls24_315,groth16,0,0
 scalar_mul_P256,bls24_317,groth16,0,0
 scalar_mul_P256,bw6_761,groth16,0,0
 scalar_mul_P256,bw6_633,groth16,0,0
-scalar_mul_P256,bn254,plonk,385060,359805
+scalar_mul_P256,bn254,plonk,355345,331788
 scalar_mul_P256,bls12_377,plonk,0,0
 scalar_mul_P256,bls12_381,plonk,0,0
 scalar_mul_P256,bls24_315,plonk,0,0
 scalar_mul_P256,bls24_317,plonk,0,0
 scalar_mul_P256,bw6_761,plonk,0,0
 scalar_mul_P256,bw6_633,plonk,0,0
-scalar_mul_secp256k1,bn254,groth16,75154,118312
+scalar_mul_secp256k1,bn254,groth16,69860,109339
 scalar_mul_secp256k1,bls12_377,groth16,0,0
 scalar_mul_secp256k1,bls12_381,groth16,0,0
 scalar_mul_secp256k1,bls24_315,groth16,0,0
 scalar_mul_secp256k1,bls24_317,groth16,0,0
 scalar_mul_secp256k1,bw6_761,groth16,0,0
 scalar_mul_secp256k1,bw6_633,groth16,0,0
-scalar_mul_secp256k1,bn254,plonk,281870,264753
+scalar_mul_secp256k1,bn254,plonk,263180,247131
 scalar_mul_secp256k1,bls12_377,plonk,0,0
 scalar_mul_secp256k1,bls12_381,plonk,0,0
 scalar_mul_secp256k1,bls24_315,plonk,0,0
diff --git a/std/internal/limbcomposition/composition.go b/std/internal/limbcomposition/composition.go
index 7fad01bf0a..b2d3211c75 100644
--- a/std/internal/limbcomposition/composition.go
+++ b/std/internal/limbcomposition/composition.go
@@ -12,9 +12,6 @@ import (
 //
 //	res = \sum_{i=0}^{len(inputs)} inputs[i] * 2^{nbBits * i}
 func Recompose(inputs []*big.Int, nbBits uint, res *big.Int) error {
-	if len(inputs) == 0 {
-		return fmt.Errorf("zero length slice input")
-	}
 	if res == nil {
 		return fmt.Errorf("result not initialized")
 	}
diff --git a/std/math/emulated/custommod.go b/std/math/emulated/custommod.go
index acc25c8eda..143ab46293 100644
--- a/std/math/emulated/custommod.go
+++ b/std/math/emulated/custommod.go
@@ -14,6 +14,10 @@ import (
 //
 // NB! circuit complexity depends on T rather on the actual length of the modulus.
 func (f *Field[T]) ModMul(a, b *Element[T], modulus *Element[T]) *Element[T] {
+	// fast path when either of the inputs is zero then result is always zero
+	if len(a.Limbs) == 0 || len(b.Limbs) == 0 {
+		return f.Zero()
+	}
 	res := f.mulMod(a, b, 0, modulus)
 	return res
 }
@@ -33,9 +37,9 @@ func (f *Field[T]) ModAdd(a, b *Element[T], modulus *Element[T]) *Element[T] {
 	for nextOverflow, err = f.addPreCond(a, b); errors.As(err, &target); nextOverflow, err = f.addPreCond(a, b) {
 		if errors.As(err, &target) {
 			if !target.reduceRight {
-				a = f.mulMod(a, f.shortOne(), 0, modulus)
+				a = f.mulMod(a, f.One(), 0, modulus)
 			} else {
-				b = f.mulMod(b, f.shortOne(), 0, modulus)
+				b = f.mulMod(b, f.One(), 0, modulus)
 			}
 		}
 	}
@@ -93,6 +97,10 @@ func (f *Field[T]) ModAssertIsEqual(a, b *Element[T], modulus *Element[T]) {
 //
 // NB! circuit complexity depends on T rather on the actual length of the modulus.
 func (f *Field[T]) ModExp(base, exp, modulus *Element[T]) *Element[T] {
+	// fasth path when the base is zero then result is always zero
+	if len(base.Limbs) == 0 {
+		return f.Zero()
+	}
 	expBts := f.ToBits(exp)
 	n := len(expBts)
 	res := f.Select(expBts[0], base, f.One())
diff --git a/std/math/emulated/element.go b/std/math/emulated/element.go
index bbad40cd2c..2603deb64c 100644
--- a/std/math/emulated/element.go
+++ b/std/math/emulated/element.go
@@ -43,14 +43,23 @@ type Element[T FieldParams] struct {
 	evaluation  frontend.Variable `gnark:"-"`
 }
 
-// ValueOf returns an Element[T] from a constant value.
-// The input is converted to *big.Int and decomposed into limbs and packed into new Element[T].
+// ValueOf returns an Element[T] from a constant value. This method is used for
+// witness assignment. For in-circuit constant assignment use the
+// [Field.NewElement] method.
+//
+// The input is converted into limbs according to the parameters of the field
+// and returned as a new [Element[T]]. Note that it returns the value, not a
+// reference, which is more convenient for witness assignment.
 func ValueOf[T FieldParams](constant interface{}) Element[T] {
+	// in this method we set the isWitness flag to true, because we do not know
+	// the width of the input value. Even though it is valid to call this method
+	// in circuit without reference to `Field`, then the canonical way would be
+	// to call [Field.NewElement] method (which would set isWitness to false).
 	if constant == nil {
-		r := newConstElement[T](0)
+		r := newConstElement[T](0, true)
 		return *r
 	}
-	r := newConstElement[T](constant)
+	r := newConstElement[T](constant, true)
 	return *r
 }
 
@@ -58,7 +67,7 @@ func ValueOf[T FieldParams](constant interface{}) Element[T] {
 // taking pointer to it. We only want to have a public method for initialising
 // an element which return a value because the user uses this only for witness
 // creation and it mess up schema parsing.
-func newConstElement[T FieldParams](v interface{}) *Element[T] {
+func newConstElement[T FieldParams](v interface{}, isWitness bool) *Element[T] {
 	var fp T
 	// convert to big.Int
 	bValue := utils.FromInterface(v)
@@ -68,9 +77,18 @@ func newConstElement[T FieldParams](v interface{}) *Element[T] {
 		bValue.Mod(&bValue, fp.Modulus())
 	}
 
-	// decompose into limbs
+	// decompose into limbs. When set with isWitness, then we do not know at
+	// compile time the width of the input, so we allocate the maximum number of
+	// limbs. However, in-circuit we already do (we set it from actual
+	// constant), thus we can allocate the exact number of limbs.
+	var nbLimbs int
+	if isWitness {
+		nbLimbs = int(fp.NbLimbs())
+	} else {
+		nbLimbs = (bValue.BitLen() + int(fp.BitsPerLimb()) - 1) / int(fp.BitsPerLimb())
+	}
 	// TODO @gbotrel use big.Int pool here
-	blimbs := make([]*big.Int, fp.NbLimbs())
+	blimbs := make([]*big.Int, nbLimbs)
 	for i := range blimbs {
 		blimbs[i] = new(big.Int)
 	}
diff --git a/std/math/emulated/element_test.go b/std/math/emulated/element_test.go
index 37389ce155..a711dda068 100644
--- a/std/math/emulated/element_test.go
+++ b/std/math/emulated/element_test.go
@@ -1425,3 +1425,75 @@ func testPolyEvalNegativeCoefficient[T FieldParams](t *testing.T) {
 	err = test.IsSolved(&PolyEvalNegativeCoefficient[T]{Inputs: make([]Element[T], nbInputs)}, assignment, testCurve.ScalarField())
 	assert.NoError(err)
 }
+
+type FastPathsCircuit[T FieldParams] struct {
+	Rand Element[T]
+	Zero Element[T]
+}
+
+func (c *FastPathsCircuit[T]) Define(api frontend.API) error {
+	f, err := NewField[T](api)
+	if err != nil {
+		return err
+	}
+	// instead of using witness values, we need to create the elements
+	// in-circuit. In witness creation we always create elements with full
+	// number of limbs.
+
+	zero := f.Zero()
+
+	// mul
+	res := f.Mul(zero, &c.Rand)
+	f.AssertIsEqual(res, &c.Zero)
+	f.AssertIsEqual(res, zero)
+	res = f.Mul(&c.Rand, zero)
+	f.AssertIsEqual(res, &c.Zero)
+	f.AssertIsEqual(res, zero)
+
+	res = f.MulMod(zero, &c.Rand)
+	f.AssertIsEqual(res, &c.Zero)
+	f.AssertIsEqual(res, zero)
+	res = f.MulMod(&c.Rand, zero)
+	f.AssertIsEqual(res, &c.Zero)
+	f.AssertIsEqual(res, zero)
+
+	res = f.MulNoReduce(zero, &c.Rand)
+	f.AssertIsEqual(res, &c.Zero)
+	f.AssertIsEqual(res, zero)
+	res = f.MulNoReduce(&c.Rand, zero)
+	f.AssertIsEqual(res, &c.Zero)
+	f.AssertIsEqual(res, zero)
+
+	// div
+	res = f.Div(zero, &c.Rand)
+	f.AssertIsEqual(res, &c.Zero)
+	f.AssertIsEqual(res, zero)
+
+	// square root
+	res = f.Sqrt(zero)
+	f.AssertIsEqual(res, &c.Zero)
+	f.AssertIsEqual(res, zero)
+
+	// exp
+	res = f.Exp(zero, &c.Rand)
+	f.AssertIsEqual(res, &c.Zero)
+	f.AssertIsEqual(res, zero)
+
+	return nil
+}
+
+func TestFasthPaths(t *testing.T) {
+	testFastPaths[Goldilocks](t)
+	testFastPaths[BN254Fr](t)
+	testFastPaths[emparams.Mod1e512](t)
+}
+
+func testFastPaths[T FieldParams](t *testing.T) {
+	assert := test.NewAssert(t)
+	var fp T
+	randVal, _ := rand.Int(rand.Reader, fp.Modulus())
+	circuit := &FastPathsCircuit[T]{}
+	assignment := &FastPathsCircuit[T]{Rand: ValueOf[T](randVal), Zero: ValueOf[T](0)}
+
+	assert.CheckCircuit(circuit, test.WithValidAssignment(assignment))
+}
diff --git a/std/math/emulated/emparams/emparams.go b/std/math/emulated/emparams/emparams.go
index 736530148c..b6ac1ce9ab 100644
--- a/std/math/emulated/emparams/emparams.go
+++ b/std/math/emulated/emparams/emparams.go
@@ -43,6 +43,11 @@ func (twelveLimbPrimeField) NbLimbs() uint     { return 12 }
 func (twelveLimbPrimeField) BitsPerLimb() uint { return 64 }
 func (twelveLimbPrimeField) IsPrime() bool     { return true }
 
+type oneLimbPrimeField struct{}
+
+func (oneLimbPrimeField) NbLimbs() uint { return 1 }
+func (oneLimbPrimeField) IsPrime() bool { return true }
+
 // Goldilocks provides type parametrization for field emulation:
 //   - limbs: 1
 //   - limb width: 64 bits
@@ -51,11 +56,9 @@ func (twelveLimbPrimeField) IsPrime() bool     { return true }
 //
 //	0xffffffff00000001 (base 16)
 //	18446744069414584321 (base 10)
-type Goldilocks struct{}
+type Goldilocks struct{ oneLimbPrimeField }
 
-func (fp Goldilocks) NbLimbs() uint     { return 1 }
 func (fp Goldilocks) BitsPerLimb() uint { return 64 }
-func (fp Goldilocks) IsPrime() bool     { return true }
 func (fp Goldilocks) Modulus() *big.Int { return goldilocks.Modulus() }
 
 // Secp256k1Fp provides type parametrization for field emulation:
@@ -366,3 +369,35 @@ func (Mod1e256) Modulus() *big.Int {
 	val, _ := new(big.Int).SetString("ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", 16)
 	return val
 }
+
+// BabyBear provides type parametrization for field emulation:
+//   - limbs: 1
+//   - limb width: 31 bits
+//
+// The prime modulus for type parametrisation is:
+//
+//	15*2^27+1
+//	0x78000001 (base 16)
+//	2013265921 (base 10)
+//
+// The field has 2-adicity of 27.
+type BabyBear struct{ oneLimbPrimeField }
+
+func (BabyBear) BitsPerLimb() uint { return 31 }
+func (BabyBear) Modulus() *big.Int { return big.NewInt(2013265921) }
+
+// KoalaBear provides type parametrization for field emulation:
+//   - limbs: 1
+//   - limb width: 31 bits
+//
+// The prime modulus for type parametrisation is:
+//
+//	2^31-2^24+1
+//	0x7f000001 (base 16)
+//	2130706433 (base 10)
+//
+// The field has 2-adicity of 24.
+type KoalaBear struct{ oneLimbPrimeField }
+
+func (KoalaBear) BitsPerLimb() uint { return 31 }
+func (KoalaBear) Modulus() *big.Int { return big.NewInt(2130706433) }
diff --git a/std/math/emulated/field.go b/std/math/emulated/field.go
index d7aae8714f..d7087974c6 100644
--- a/std/math/emulated/field.go
+++ b/std/math/emulated/field.go
@@ -129,14 +129,17 @@ func (f *Field[T]) NewElement(v interface{}) *Element[T] {
 	if e, ok := v.([]frontend.Variable); ok {
 		return f.packLimbs(e, true)
 	}
-	c := ValueOf[T](v)
-	return &c
+	// the input was not a variable, so it must be a constant. Create a new
+	// element from it while setting isWitness flag to false. This ensures that
+	// we use the minimal number of limbs necessary.
+	c := newConstElement[T](v, false)
+	return c
 }
 
 // Zero returns zero as a constant.
 func (f *Field[T]) Zero() *Element[T] {
 	f.zeroConstOnce.Do(func() {
-		f.zeroConst = newConstElement[T](0)
+		f.zeroConst = f.newInternalElement([]frontend.Variable{}, 0)
 	})
 	return f.zeroConst
 }
@@ -144,23 +147,15 @@ func (f *Field[T]) Zero() *Element[T] {
 // One returns one as a constant.
 func (f *Field[T]) One() *Element[T] {
 	f.oneConstOnce.Do(func() {
-		f.oneConst = newConstElement[T](1)
+		f.oneConst = f.newInternalElement([]frontend.Variable{1}, 0)
 	})
 	return f.oneConst
 }
 
-// shortOne returns one as a constant stored in a single limb.
-func (f *Field[T]) shortOne() *Element[T] {
-	f.shortOneConstOnce.Do(func() {
-		f.shortOneConst = f.newInternalElement([]frontend.Variable{1}, 0)
-	})
-	return f.shortOneConst
-}
-
 // Modulus returns the modulus of the emulated ring as a constant.
 func (f *Field[T]) Modulus() *Element[T] {
 	f.nConstOnce.Do(func() {
-		f.nConst = newConstElement[T](f.fParams.Modulus())
+		f.nConst = newConstElement[T](f.fParams.Modulus(), false)
 	})
 	return f.nConst
 }
@@ -168,7 +163,7 @@ func (f *Field[T]) Modulus() *Element[T] {
 // modulusPrev returns modulus-1 as a constant.
 func (f *Field[T]) modulusPrev() *Element[T] {
 	f.nprevConstOnce.Do(func() {
-		f.nprevConst = newConstElement[T](new(big.Int).Sub(f.fParams.Modulus(), big.NewInt(1)))
+		f.nprevConst = newConstElement[T](new(big.Int).Sub(f.fParams.Modulus(), big.NewInt(1)), false)
 	})
 	return f.nprevConst
 }
diff --git a/std/math/emulated/field_mul.go b/std/math/emulated/field_mul.go
index d6e3bfdfb0..224ffb4cc9 100644
--- a/std/math/emulated/field_mul.go
+++ b/std/math/emulated/field_mul.go
@@ -214,11 +214,15 @@ func (f *Field[T]) mulMod(a, b *Element[T], _ uint, p *Element[T]) *Element[T] {
 
 // checkZero creates multiplication check a * 1 = 0 + k*p.
 func (f *Field[T]) checkZero(a *Element[T], p *Element[T]) {
+	// fast path - the result is on zero limbs. This means that it is constant zero
+	if len(a.Limbs) == 0 {
+		return
+	}
 	// the method works similarly to mulMod, but we know that we are multiplying
 	// by one and expected result should be zero.
 	f.enforceWidthConditional(a)
 	f.enforceWidthConditional(p)
-	b := f.shortOne()
+	b := f.One()
 	k, r, c, err := f.callMulHint(a, b, false, p)
 	if err != nil {
 		panic(err)
@@ -454,6 +458,10 @@ func mulHint(field *big.Int, inputs, outputs []*big.Int) error {
 // For multiplying by a constant, use [Field[T].MulConst] method which is more
 // efficient.
 func (f *Field[T]) Mul(a, b *Element[T]) *Element[T] {
+	// fast path - if one of the inputs is on zero limbs (it is zero), then the result is also zero
+	if len(a.Limbs) == 0 || len(b.Limbs) == 0 {
+		return f.Zero()
+	}
 	return f.reduceAndOp(func(a, b *Element[T], u uint) *Element[T] { return f.mulMod(a, b, u, nil) }, f.mulPreCond, a, b)
 }
 
@@ -462,6 +470,10 @@ func (f *Field[T]) Mul(a, b *Element[T]) *Element[T] {
 //
 // Equivalent to [Field[T].Mul], kept for backwards compatibility.
 func (f *Field[T]) MulMod(a, b *Element[T]) *Element[T] {
+	// fast path - if one of the inputs is on zero limbs (it is zero), then the result is also zero
+	if len(a.Limbs) == 0 || len(b.Limbs) == 0 {
+		return f.Zero()
+	}
 	return f.reduceAndOp(func(a, b *Element[T], u uint) *Element[T] { return f.mulMod(a, b, u, nil) }, f.mulPreCond, a, b)
 }
 
@@ -471,6 +483,9 @@ func (f *Field[T]) MulMod(a, b *Element[T]) *Element[T] {
 // general [Field[T].Mul] or [Field[T].MulMod] with creating new Element from
 // the constant on-the-fly.
 func (f *Field[T]) MulConst(a *Element[T], c *big.Int) *Element[T] {
+	if len(a.Limbs) == 0 {
+		return f.Zero()
+	}
 	switch c.Sign() {
 	case -1:
 		f.MulConst(f.Neg(a), new(big.Int).Neg(c))
@@ -485,7 +500,7 @@ func (f *Field[T]) MulConst(a *Element[T], c *big.Int) *Element[T] {
 		func(a, _ *Element[T], u uint) *Element[T] {
 			if ba, aConst := f.constantValue(a); aConst {
 				ba.Mul(ba, c)
-				return newConstElement[T](ba)
+				return newConstElement[T](ba, false)
 			}
 			limbs := make([]frontend.Variable, len(a.Limbs))
 			for i := range a.Limbs {
@@ -522,6 +537,10 @@ func (f *Field[T]) mulPreCond(a, b *Element[T]) (nextOverflow uint, err error) {
 // the field order. The number of limbs of the returned element depends on the
 // number of limbs of the inputs.
 func (f *Field[T]) MulNoReduce(a, b *Element[T]) *Element[T] {
+	// fast path - if one of the inputs is on zero limbs (it is zero), then the result is also zero
+	if len(a.Limbs) == 0 || len(b.Limbs) == 0 {
+		return f.Zero()
+	}
 	return f.reduceAndOp(f.mulNoReduce, f.mulPreCond, a, b)
 }
 
@@ -541,6 +560,10 @@ func (f *Field[T]) mulNoReduce(a, b *Element[T], nextoverflow uint) *Element[T]
 // Exp computes base^exp modulo the field order. The returned Element has default
 // number of limbs and zero overflow.
 func (f *Field[T]) Exp(base, exp *Element[T]) *Element[T] {
+	// fast path - if the base is zero, then the result is also zero
+	if len(base.Limbs) == 0 {
+		return f.Zero()
+	}
 	expBts := f.ToBits(exp)
 	n := len(expBts)
 	res := f.Select(expBts[0], base, f.One())
diff --git a/std/math/emulated/field_ops.go b/std/math/emulated/field_ops.go
index 9ee181b7fb..06b52c541a 100644
--- a/std/math/emulated/field_ops.go
+++ b/std/math/emulated/field_ops.go
@@ -11,6 +11,10 @@ import (
 
 // Div computes a/b and returns it. It uses [DivHint] as a hint function.
 func (f *Field[T]) Div(a, b *Element[T]) *Element[T] {
+	// fast path when dividing by 0
+	if len(a.Limbs) == 0 {
+		return f.Zero()
+	}
 	return f.reduceAndOp(f.div, f.divPreCond, a, b)
 }
 
@@ -70,6 +74,10 @@ func (f *Field[T]) inverse(a, _ *Element[T], _ uint) *Element[T] {
 
 // Sqrt computes square root of a and returns it. It uses [SqrtHint].
 func (f *Field[T]) Sqrt(a *Element[T]) *Element[T] {
+	// fast path when input is zero
+	if len(a.Limbs) == 0 {
+		return f.Zero()
+	}
 	return f.reduceAndOp(f.sqrt, f.sqrtPreCond, a, nil)
 }
 
@@ -116,7 +124,7 @@ func (f *Field[T]) add(a, b *Element[T], nextOverflow uint) *Element[T] {
 	bb, bConst := f.constantValue(b)
 	if aConst && bConst {
 		ba.Add(ba, bb).Mod(ba, f.fParams.Modulus())
-		return newConstElement[T](ba)
+		return newConstElement[T](ba, false)
 	}
 
 	nbLimbs := max(len(a.Limbs), len(b.Limbs))
@@ -184,13 +192,13 @@ func (f *Field[T]) sub(a, b *Element[T], nextOverflow uint) *Element[T] {
 	bb, bConst := f.constantValue(b)
 	if aConst && bConst {
 		ba.Sub(ba, bb).Mod(ba, f.fParams.Modulus())
-		return newConstElement[T](ba)
+		return newConstElement[T](ba, false)
 	}
 
 	// first we have to compute padding to ensure that the subtraction does not
 	// underflow.
 	var fp T
-	nbLimbs := max(len(a.Limbs), len(b.Limbs))
+	nbLimbs := max(len(a.Limbs), len(b.Limbs), int(fp.NbLimbs()))
 	limbs := make([]frontend.Variable, nbLimbs)
 	padLimbs := subPadding(fp.Modulus(), fp.BitsPerLimb(), b.overflow, uint(nbLimbs))
 	for i := range limbs {
diff --git a/std/math/emulated/params.go b/std/math/emulated/params.go
index 2b0dc9d179..0b12de104c 100644
--- a/std/math/emulated/params.go
+++ b/std/math/emulated/params.go
@@ -18,6 +18,7 @@ import (
 //   - [P256Fp] and [P256Fr]
 //   - [P384Fp] and [P384Fr]
 //   - [STARKCurveFp] and [STARKCurveFr]
+//   - [BabyBear] and [KoalaBear]
 type FieldParams interface {
 	NbLimbs() uint     // number of limbs to represent field element
 	BitsPerLimb() uint // number of bits per limb. Top limb may contain less than limbSize bits.
@@ -42,4 +43,6 @@ type (
 	BW6761Fr     = emparams.BW6761Fr
 	STARKCurveFp = emparams.STARKCurveFp
 	STARKCurveFr = emparams.STARKCurveFr
+	BabyBear     = emparams.BabyBear
+	KoalaBear    = emparams.KoalaBear
 )
diff --git a/std/math/emulated/subtraction_padding.go b/std/math/emulated/subtraction_padding.go
index e30298e9aa..ae4ff8ddc6 100644
--- a/std/math/emulated/subtraction_padding.go
+++ b/std/math/emulated/subtraction_padding.go
@@ -26,6 +26,12 @@ func subPadding(modulus *big.Int, bitsPerLimbs uint, overflow uint, nbLimbs uint
 	if modulus.Cmp(big.NewInt(0)) == 0 {
 		panic("modulus is zero")
 	}
+	// we ensure that the number of padding limbs is sufficient to represent the
+	// modulus
+	requiredLimbs := (uint(modulus.BitLen()) + bitsPerLimbs - 1) / bitsPerLimbs
+	if nbLimbs < requiredLimbs {
+		nbLimbs = requiredLimbs
+	}
 	// first, we build a number nLimbs, such that nLimbs > b;
 	// here b is defined by its bounds, that is b is an element with nbLimbs of (bitsPerLimbs+overflow)
 	// so a number nLimbs > b, is simply taking the next power of 2 over this bound .