diff --git a/field/babybear/extensions/doc.go b/field/babybear/extensions/doc.go
new file mode 100644
index 0000000000..7557b8a6fd
--- /dev/null
+++ b/field/babybear/extensions/doc.go
@@ -0,0 +1,6 @@
+// Package extensions implements the fields arithmetic of the 𝔽r² and 𝔽r⁴
+// extensions of the Koala-bear field.
+//
+//	𝔽r²[u] = 𝔽r/u²-11
+//	𝔽r⁴[v] = 𝔽r²/v²-u
+package extensions
diff --git a/field/babybear/extensions/e2.go b/field/babybear/extensions/e2.go
new file mode 100644
index 0000000000..b602ba907e
--- /dev/null
+++ b/field/babybear/extensions/e2.go
@@ -0,0 +1,278 @@
+package extensions
+
+import (
+	"math/big"
+
+	fr "github.com/consensys/gnark-crypto/field/babybear"
+)
+
+// E2 is a degree two finite field extension of fr.Element
+type E2 struct {
+	A0, A1 fr.Element
+}
+
+// Equal returns true if z equals x, false otherwise
+func (z *E2) Equal(x *E2) bool {
+	return z.A0.Equal(&x.A0) && z.A1.Equal(&x.A1)
+}
+
+// Cmp compares (lexicographic order) z and x and returns:
+//
+//	-1 if z <  x
+//	 0 if z == x
+//	+1 if z >  x
+func (z *E2) Cmp(x *E2) int {
+	if a1 := z.A1.Cmp(&x.A1); a1 != 0 {
+		return a1
+	}
+	return z.A0.Cmp(&x.A0)
+}
+
+// LexicographicallyLargest returns true if this element is strictly lexicographically
+// larger than its negation, false otherwise
+func (z *E2) LexicographicallyLargest() bool {
+	// adapted from github.com/zkcrypto/bls12_381
+	if z.A1.IsZero() {
+		return z.A0.LexicographicallyLargest()
+	}
+	return z.A1.LexicographicallyLargest()
+}
+
+// SetString sets a E2 element from strings
+func (z *E2) SetString(s1, s2 string) *E2 {
+	z.A0.SetString(s1)
+	z.A1.SetString(s2)
+	return z
+}
+
+// SetZero sets an E2 elmt to zero
+func (z *E2) SetZero() *E2 {
+	z.A0.SetZero()
+	z.A1.SetZero()
+	return z
+}
+
+// Set sets an E2 from x
+func (z *E2) Set(x *E2) *E2 {
+	z.A0 = x.A0
+	z.A1 = x.A1
+	return z
+}
+
+// SetOne sets z to 1 in Montgomery form and returns z
+func (z *E2) SetOne() *E2 {
+	z.A0.SetOne()
+	z.A1.SetZero()
+	return z
+}
+
+// SetRandom sets a0 and a1 to random values
+func (z *E2) SetRandom() (*E2, error) {
+	if _, err := z.A0.SetRandom(); err != nil {
+		return nil, err
+	}
+	if _, err := z.A1.SetRandom(); err != nil {
+		return nil, err
+	}
+	return z, nil
+}
+
+// IsZero returns true if z is zero, false otherwise
+func (z *E2) IsZero() bool {
+	return z.A0.IsZero() && z.A1.IsZero()
+}
+
+// IsOne returns true if z is one, false otherwise
+func (z *E2) IsOne() bool {
+	return z.A0.IsOne() && z.A1.IsZero()
+}
+
+// Add adds two elements of E2
+func (z *E2) Add(x, y *E2) *E2 {
+	z.A0.Add(&x.A0, &y.A0)
+	z.A1.Add(&x.A1, &y.A1)
+	return z
+}
+
+// Sub subtracts two elements of E2
+func (z *E2) Sub(x, y *E2) *E2 {
+	z.A0.Sub(&x.A0, &y.A0)
+	z.A1.Sub(&x.A1, &y.A1)
+	return z
+}
+
+// Double doubles an E2 element
+func (z *E2) Double(x *E2) *E2 {
+	z.A0.Double(&x.A0)
+	z.A1.Double(&x.A1)
+	return z
+}
+
+// Neg negates an E2 element
+func (z *E2) Neg(x *E2) *E2 {
+	z.A0.Neg(&x.A0)
+	z.A1.Neg(&x.A1)
+	return z
+}
+
+// String implements Stringer interface for fancy printing
+func (z *E2) String() string {
+	return z.A0.String() + "+" + z.A1.String() + "*u"
+}
+
+// MulByElement multiplies an element in E2 by an element in fr
+func (z *E2) MulByElement(x *E2, y *fr.Element) *E2 {
+	var yCopy fr.Element
+	yCopy.Set(y)
+	z.A0.Mul(&x.A0, &yCopy)
+	z.A1.Mul(&x.A1, &yCopy)
+	return z
+}
+
+// Conjugate conjugates an element in E2
+func (z *E2) Conjugate(x *E2) *E2 {
+	z.A0 = x.A0
+	z.A1.Neg(&x.A1)
+	return z
+}
+
+// Halve sets z to z / 2
+func (z *E2) Halve() {
+	z.A0.Halve()
+	z.A1.Halve()
+}
+
+// Legendre returns the Legendre symbol of z
+func (z *E2) Legendre() int {
+	var n fr.Element
+	z.norm(&n)
+	return n.Legendre()
+}
+
+// Exp sets z=xᵏ (mod q²) and returns it
+func (z *E2) Exp(x E2, k *big.Int) *E2 {
+	if k.IsUint64() && k.Uint64() == 0 {
+		return z.SetOne()
+	}
+
+	e := k
+	if k.Sign() == -1 {
+		// negative k, we invert
+		// if k < 0: xᵏ (mod q²) == (x⁻¹)ᵏ (mod q²)
+		x.Inverse(&x)
+
+		// we negate k in a temp big.Int since
+		// Int.Bit(_) of k and -k is different
+		e = bigIntPool.Get().(*big.Int)
+		defer bigIntPool.Put(e)
+		e.Neg(k)
+	}
+
+	z.SetOne()
+	b := e.Bytes()
+	for i := 0; i < len(b); i++ {
+		w := b[i]
+		for j := 0; j < 8; j++ {
+			z.Square(z)
+			if (w & (0b10000000 >> j)) != 0 {
+				z.Mul(z, &x)
+			}
+		}
+	}
+
+	return z
+}
+
+// Sqrt sets z to the square root of and returns z
+// The function does not test whether the square root
+// exists or not, it's up to the caller to call
+// Legendre beforehand.
+// cf https://eprint.iacr.org/2012/685.pdf (algo 10)
+func (z *E2) Sqrt(x *E2) *E2 {
+
+	// precomputation
+	var b, c, d, e, f, x0 E2
+	var _b, o fr.Element
+
+	// c must be a non square (p = 1 mod 4)
+	c.A1.SetOne()
+
+	q := fr.Modulus()
+	var exp, one big.Int
+	one.SetUint64(1)
+	exp.Set(q).Sub(&exp, &one).Rsh(&exp, 1)
+	d.Exp(c, &exp)
+	e.Mul(&d, &c).Inverse(&e)
+	f.Mul(&d, &c).Square(&f)
+
+	// computation
+	exp.Rsh(&exp, 1)
+	b.Exp(*x, &exp)
+	b.norm(&_b)
+	o.SetOne()
+	if _b.Equal(&o) {
+		x0.Square(&b).Mul(&x0, x)
+		_b.Set(&x0.A0).Sqrt(&_b)
+		z.Conjugate(&b).MulByElement(z, &_b)
+		return z
+	}
+	x0.Square(&b).Mul(&x0, x).Mul(&x0, &f)
+	_b.Set(&x0.A0).Sqrt(&_b)
+	z.Conjugate(&b).MulByElement(z, &_b).Mul(z, &e)
+
+	return z
+}
+
+// BatchInvertE2 returns a new slice with every element in a inverted.
+// It uses Montgomery batch inversion trick.
+//
+// if a[i] == 0, returns result[i] = a[i]
+func BatchInvertE2(a []E2) []E2 {
+	res := make([]E2, len(a))
+	if len(a) == 0 {
+		return res
+	}
+
+	zeroes := make([]bool, len(a))
+	var accumulator E2
+	accumulator.SetOne()
+
+	for i := 0; i < len(a); i++ {
+		if a[i].IsZero() {
+			zeroes[i] = true
+			continue
+		}
+		res[i].Set(&accumulator)
+		accumulator.Mul(&accumulator, &a[i])
+	}
+
+	accumulator.Inverse(&accumulator)
+
+	for i := len(a) - 1; i >= 0; i-- {
+		if zeroes[i] {
+			continue
+		}
+		res[i].Mul(&res[i], &accumulator)
+		accumulator.Mul(&accumulator, &a[i])
+	}
+
+	return res
+}
+
+// Select is conditional move.
+// If cond = 0, it sets z to caseZ and returns it. otherwise caseNz.
+func (z *E2) Select(cond int, caseZ *E2, caseNz *E2) *E2 {
+	//Might be able to save a nanosecond or two by an aggregate implementation
+
+	z.A0.Select(cond, &caseZ.A0, &caseNz.A0)
+	z.A1.Select(cond, &caseZ.A1, &caseNz.A1)
+
+	return z
+}
+
+// Div divides an element in E2 by an element in E2
+func (z *E2) Div(x *E2, y *E2) *E2 {
+	var r E2
+	r.Inverse(y).Mul(x, &r)
+	return z.Set(&r)
+}
diff --git a/field/babybear/extensions/e2_babybear.go b/field/babybear/extensions/e2_babybear.go
new file mode 100644
index 0000000000..dd45cd43b9
--- /dev/null
+++ b/field/babybear/extensions/e2_babybear.go
@@ -0,0 +1,87 @@
+package extensions
+
+import (
+	fr "github.com/consensys/gnark-crypto/field/babybear"
+)
+
+// MulBy11 x *= 11 (mod q)
+func MulBy11(x *fr.Element) {
+	var y fr.Element
+	y.SetUint64(11)
+	x.Mul(x, &y)
+
+}
+
+// Mul sets z to the E2-product of x,y, returns z
+func (z *E2) Mul(x, y *E2) *E2 {
+	var a, b, c fr.Element
+	a.Add(&x.A0, &x.A1)
+	b.Add(&y.A0, &y.A1)
+	a.Mul(&a, &b)
+	b.Mul(&x.A0, &y.A0)
+	c.Mul(&x.A1, &y.A1)
+	z.A1.Sub(&a, &b).Sub(&z.A1, &c)
+	MulBy11(&c)
+	z.A0.Add(&b, &c)
+	return z
+}
+
+// Square sets z to the E2-product of x,x returns z
+func (z *E2) Square(x *E2) *E2 {
+	var a, b, c fr.Element
+	a.Mul(&x.A0, &x.A1).Double(&a)
+	c.Square(&x.A0)
+	b.Square(&x.A1)
+	MulBy11(&b)
+	z.A0.Add(&c, &b)
+	z.A1 = a
+	return z
+}
+
+// MulByNonResidue multiplies a E2 by (0,1)
+func (z *E2) MulByNonResidue(x *E2) *E2 {
+	a := x.A0
+	b := x.A1 // fetching x.A1 in the function below is slower
+	MulBy11(&b)
+	z.A0 = b
+	z.A1 = a
+	return z
+}
+
+// MulByNonResidueInv multiplies a E2 by (0,1)^{-1}
+func (z *E2) MulByNonResidueInv(x *E2) *E2 {
+	a := x.A1
+	// 1/11 mod r
+	var elevenInv fr.Element
+	elevenInv.SetUint64(549072524)
+	z.A1.Mul(&x.A0, &elevenInv)
+	z.A0 = a
+	return z
+}
+
+// Inverse sets z to the E2-inverse of x, returns z
+func (z *E2) Inverse(x *E2) *E2 {
+	// Algorithm 8 from https://eprint.iacr.org/2010/354.pdf
+	var t0, t1, tmp fr.Element
+	a := &x.A0 // creating the buffers a, b is faster than querying &x.A0, &x.A1 in the functions call below
+	b := &x.A1
+	t0.Square(a)
+	t1.Square(b)
+	tmp.Set(&t1)
+	MulBy11(&tmp)
+	t0.Sub(&t0, &tmp)
+	t1.Inverse(&t0)
+	z.A0.Mul(a, &t1)
+	z.A1.Mul(b, &t1).Neg(&z.A1)
+
+	return z
+}
+
+// norm sets x to the norm of z
+func (z *E2) norm(x *fr.Element) {
+	var tmp fr.Element
+	x.Square(&z.A1)
+	tmp.Set(x)
+	MulBy11(&tmp)
+	x.Square(&z.A0).Sub(x, &tmp)
+}
diff --git a/field/babybear/extensions/e2_test.go b/field/babybear/extensions/e2_test.go
new file mode 100644
index 0000000000..103fec029c
--- /dev/null
+++ b/field/babybear/extensions/e2_test.go
@@ -0,0 +1,510 @@
+package extensions
+
+import (
+	"crypto/rand"
+	"testing"
+
+	fr "github.com/consensys/gnark-crypto/field/babybear"
+	"github.com/leanovate/gopter"
+	"github.com/leanovate/gopter/prop"
+)
+
+// ------------------------------------------------------------
+// tests
+
+const (
+	nbFuzzShort = 10
+	nbFuzz      = 50
+)
+
+func TestE2ReceiverIsOperand(t *testing.T) {
+
+	t.Parallel()
+	parameters := gopter.DefaultTestParameters()
+	if testing.Short() {
+		parameters.MinSuccessfulTests = nbFuzzShort
+	} else {
+		parameters.MinSuccessfulTests = nbFuzz
+	}
+
+	properties := gopter.NewProperties(parameters)
+
+	genA := GenE2()
+	genB := GenE2()
+	genfr := GenFr()
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (addition) should output the same result", prop.ForAll(
+		func(a, b *E2) bool {
+			var c, d E2
+			d.Set(a)
+			c.Add(a, b)
+			a.Add(a, b)
+			b.Add(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (sub) should output the same result", prop.ForAll(
+		func(a, b *E2) bool {
+			var c, d E2
+			d.Set(a)
+			c.Sub(a, b)
+			a.Sub(a, b)
+			b.Sub(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (mul) should output the same result", prop.ForAll(
+		func(a, b *E2) bool {
+			var c, d E2
+			d.Set(a)
+			c.Mul(a, b)
+			a.Mul(a, b)
+			b.Mul(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (square) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Square(a)
+			a.Square(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (neg) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Neg(a)
+			a.Neg(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (double) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Double(a)
+			a.Double(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (mul by non residue) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.MulByNonResidue(a)
+			a.MulByNonResidue(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (mul by non residue inverse) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.MulByNonResidueInv(a)
+			a.MulByNonResidueInv(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (Inverse) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Inverse(a)
+			a.Inverse(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (Conjugate) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Conjugate(a)
+			a.Conjugate(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (mul by element) should output the same result", prop.ForAll(
+		func(a *E2, b fr.Element) bool {
+			var c E2
+			c.MulByElement(a, &b)
+			a.MulByElement(a, &b)
+			return a.Equal(&c)
+		},
+		genA,
+		genfr,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (Sqrt) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b, c, d, s E2
+
+			s.Square(a)
+			a.Set(&s)
+			b.Set(&s)
+
+			a.Sqrt(a)
+			b.Sqrt(&b)
+
+			c.Square(a)
+			d.Square(&b)
+			return c.Equal(&d)
+		},
+		genA,
+	))
+
+	properties.TestingRun(t, gopter.ConsoleReporter(false))
+
+}
+
+func TestE2MulMaxed(t *testing.T) {
+	// let's pick a and b, with maxed A0 and A1
+	var a, b E2
+	frMaxValue := fr.Element{
+		2013265921,
+	}
+	frMaxValue[0]--
+
+	a.A0 = frMaxValue
+	a.A1 = frMaxValue
+	b.A0 = frMaxValue
+	b.A1 = frMaxValue
+
+	var c, d E2
+	d.Inverse(&b)
+	c.Set(&a)
+	c.Mul(&c, &b).Mul(&c, &d)
+	if !c.Equal(&a) {
+		t.Fatal("mul with max fr failed")
+	}
+}
+
+func TestE2Ops(t *testing.T) {
+
+	t.Parallel()
+	parameters := gopter.DefaultTestParameters()
+	if testing.Short() {
+		parameters.MinSuccessfulTests = nbFuzzShort
+	} else {
+		parameters.MinSuccessfulTests = nbFuzz
+	}
+
+	properties := gopter.NewProperties(parameters)
+
+	genA := GenE2()
+	genB := GenE2()
+	genfr := GenFr()
+
+	properties.Property("[BABYBEAR] sub & add should leave an element invariant", prop.ForAll(
+		func(a, b *E2) bool {
+			var c E2
+			c.Set(a)
+			c.Add(&c, b).Sub(&c, b)
+			return c.Equal(a)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[BABYBEAR] mul & inverse should leave an element invariant", prop.ForAll(
+		func(a, b *E2) bool {
+			var c, d E2
+			d.Inverse(b)
+			c.Set(a)
+			c.Mul(&c, b).Mul(&c, &d)
+			return c.Equal(a)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[BABYBEAR] BatchInvertE2 should output the same result as Inverse", prop.ForAll(
+		func(a, b, c *E2) bool {
+
+			batch := BatchInvertE2([]E2{*a, *b, *c})
+			a.Inverse(a)
+			b.Inverse(b)
+			c.Inverse(c)
+			return a.Equal(&batch[0]) && b.Equal(&batch[1]) && c.Equal(&batch[2])
+		},
+		genA,
+		genA,
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] inverse twice should leave an element invariant", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Inverse(a).Inverse(&b)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] neg twice should leave an element invariant", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Neg(a).Neg(&b)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] square and mul should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b, c E2
+			b.Mul(a, a)
+			c.Square(a)
+			return b.Equal(&c)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] MulByElement MulByElement inverse should leave an element invariant", prop.ForAll(
+		func(a *E2, b fr.Element) bool {
+			var c E2
+			var d fr.Element
+			d.Inverse(&b)
+			c.MulByElement(a, &b).MulByElement(&c, &d)
+			return c.Equal(a)
+		},
+		genA,
+		genfr,
+	))
+
+	properties.Property("[BABYBEAR] Double and mul by 2 should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			var c fr.Element
+			c.SetUint64(2)
+			b.Double(a)
+			a.MulByElement(a, &c)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Mulbynonres mulbynonresinv should leave the element invariant", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.MulByNonResidue(a).MulByNonResidueInv(&b)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] a + pi(a), a-pi(a) should be real", prop.ForAll(
+		func(a *E2) bool {
+			var b, c, d E2
+			var e, f fr.Element
+			b.Conjugate(a)
+			c.Add(a, &b)
+			d.Sub(a, &b)
+			e.Double(&a.A0)
+			f.Double(&a.A1)
+			return c.A1.IsZero() && d.A0.IsZero() && e.Equal(&c.A0) && f.Equal(&d.A1)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Legendre on square should output 1", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Square(a)
+			c := b.Legendre()
+			return c == 1
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] square(sqrt) should leave an element invariant", prop.ForAll(
+		func(a *E2) bool {
+			var b, c, d, e E2
+			b.Square(a)
+			c.Sqrt(&b)
+			d.Square(&c)
+			e.Neg(a)
+			return (c.Equal(a) || c.Equal(&e)) && d.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] neg(E2) == neg(E2.A0, E2.A1)", prop.ForAll(
+		func(a *E2) bool {
+			var b, c E2
+			b.Neg(a)
+			c.A0.Neg(&a.A0)
+			c.A1.Neg(&a.A1)
+			return c.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Cmp and LexicographicallyLargest should be consistent", prop.ForAll(
+		func(a *E2) bool {
+			var negA E2
+			negA.Neg(a)
+			cmpResult := a.Cmp(&negA)
+			lResult := a.LexicographicallyLargest()
+			if lResult && cmpResult == 1 {
+				return true
+			}
+			if !lResult && cmpResult != 1 {
+				return true
+			}
+			return false
+		},
+		genA,
+	))
+
+	properties.TestingRun(t, gopter.ConsoleReporter(false))
+
+}
+
+// ------------------------------------------------------------
+// benches
+
+func BenchmarkE2Add(b *testing.B) {
+	var a, c E2
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Add(&a, &c)
+	}
+}
+
+func BenchmarkE2Sub(b *testing.B) {
+	var a, c E2
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Sub(&a, &c)
+	}
+}
+
+func BenchmarkE2Mul(b *testing.B) {
+	var a, c E2
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Mul(&a, &c)
+	}
+}
+
+func BenchmarkE2MulByElement(b *testing.B) {
+	var a E2
+	var c fr.Element
+	_, _ = c.SetRandom()
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.MulByElement(&a, &c)
+	}
+}
+
+func BenchmarkE2Square(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Square(&a)
+	}
+}
+
+func BenchmarkE2Sqrt(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Sqrt(&a)
+	}
+}
+
+func BenchmarkE2Exp(b *testing.B) {
+	var x E2
+	_, _ = x.SetRandom()
+	b1, _ := rand.Int(rand.Reader, fr.Modulus())
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		x.Exp(x, b1)
+	}
+}
+
+func BenchmarkE2Inverse(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Inverse(&a)
+	}
+}
+
+func BenchmarkE2MulNonRes(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.MulByNonResidue(&a)
+	}
+}
+
+func BenchmarkE2MulNonResInv(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.MulByNonResidueInv(&a)
+	}
+}
+
+func BenchmarkE2Conjugate(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Conjugate(&a)
+	}
+}
+
+func TestE2Div(t *testing.T) {
+
+	parameters := gopter.DefaultTestParameters()
+	properties := gopter.NewProperties(parameters)
+
+	genA := GenE2()
+	genB := GenE2()
+
+	properties.Property("[BABYBEAR] dividing then multiplying by the same element does nothing", prop.ForAll(
+		func(a, b *E2) bool {
+			var c E2
+			c.Div(a, b)
+			c.Mul(&c, b)
+			return c.Equal(a)
+		},
+		genA,
+		genB,
+	))
+
+	properties.TestingRun(t, gopter.ConsoleReporter(false))
+}
diff --git a/field/babybear/extensions/e4.go b/field/babybear/extensions/e4.go
new file mode 100644
index 0000000000..d4af517b14
--- /dev/null
+++ b/field/babybear/extensions/e4.go
@@ -0,0 +1,341 @@
+package extensions
+
+import (
+	"math/big"
+
+	fr "github.com/consensys/gnark-crypto/field/babybear"
+)
+
+// E4 is a degree two finite field extension of fr2
+type E4 struct {
+	B0, B1 E2
+}
+
+// Equal returns true if z equals x, false otherwise
+func (z *E4) Equal(x *E4) bool {
+	return z.B0.Equal(&x.B0) && z.B1.Equal(&x.B1)
+}
+
+// Cmp compares (lexicographic order) z and x and returns:
+//
+//	-1 if z <  x
+//	 0 if z == x
+//	+1 if z >  x
+func (z *E4) Cmp(x *E4) int {
+	if a1 := z.B1.Cmp(&x.B1); a1 != 0 {
+		return a1
+	}
+	return z.B0.Cmp(&x.B0)
+}
+
+// LexicographicallyLargest returns true if this element is strictly lexicographically
+// larger than its negation, false otherwise
+func (z *E4) LexicographicallyLargest() bool {
+	// adapted from github.com/zkcrypto/bls12_381
+	if z.B1.IsZero() {
+		return z.B0.LexicographicallyLargest()
+	}
+	return z.B1.LexicographicallyLargest()
+}
+
+// String puts E4 in string form
+func (z *E4) String() string {
+	return (z.B0.String() + "+(" + z.B1.String() + ")*v")
+}
+
+// SetString sets a E4 from string
+func (z *E4) SetString(s0, s1, s2, s3 string) *E4 {
+	z.B0.SetString(s0, s1)
+	z.B1.SetString(s2, s3)
+	return z
+}
+
+// Set copies x into z and returns z
+func (z *E4) Set(x *E4) *E4 {
+	z.B0 = x.B0
+	z.B1 = x.B1
+	return z
+}
+
+// SetZero sets an E4 elmt to zero
+func (z *E4) SetZero() *E4 {
+	z.B0.SetZero()
+	z.B1.SetZero()
+	return z
+}
+
+// SetOne sets z to 1 in Montgomery form and returns z
+func (z *E4) SetOne() *E4 {
+	*z = E4{}
+	z.B0.A0.SetOne()
+	return z
+}
+
+// MulByElement multiplies an element in E4 by an element in fr
+func (z *E4) MulByElement(x *E4, y *fr.Element) *E4 {
+	var yCopy fr.Element
+	yCopy.Set(y)
+	z.B0.MulByElement(&x.B0, &yCopy)
+	z.B1.MulByElement(&x.B1, &yCopy)
+	return z
+}
+
+// MulByE2 multiplies an element in E4 by an element in E2
+func (z *E4) MulByE2(x *E4, y *E2) *E4 {
+	var yCopy E2
+	yCopy.Set(y)
+	z.B0.Mul(&x.B0, &yCopy)
+	z.B1.Mul(&x.B1, &yCopy)
+	return z
+}
+
+// Add sets z=x+y in E4 and returns z
+func (z *E4) Add(x, y *E4) *E4 {
+	z.B0.Add(&x.B0, &y.B0)
+	z.B1.Add(&x.B1, &y.B1)
+	return z
+}
+
+// Sub sets z to x-y and returns z
+func (z *E4) Sub(x, y *E4) *E4 {
+	z.B0.Sub(&x.B0, &y.B0)
+	z.B1.Sub(&x.B1, &y.B1)
+	return z
+}
+
+// Double sets z=2*x and returns z
+func (z *E4) Double(x *E4) *E4 {
+	z.B0.Double(&x.B0)
+	z.B1.Double(&x.B1)
+	return z
+}
+
+// Neg negates an E4 element
+func (z *E4) Neg(x *E4) *E4 {
+	z.B0.Neg(&x.B0)
+	z.B1.Neg(&x.B1)
+	return z
+}
+
+// SetRandom used only in tests
+func (z *E4) SetRandom() (*E4, error) {
+	if _, err := z.B0.SetRandom(); err != nil {
+		return nil, err
+	}
+	if _, err := z.B1.SetRandom(); err != nil {
+		return nil, err
+	}
+	return z, nil
+}
+
+// IsZero returns true if z is zero, false otherwise
+func (z *E4) IsZero() bool {
+	return z.B0.IsZero() && z.B1.IsZero()
+}
+
+// IsOne returns true if z is one, false otherwise
+func (z *E4) IsOne() bool {
+	return z.B0.IsOne() && z.B1.IsZero()
+}
+
+// MulByNonResidue mul x by (0,1)
+func (z *E4) MulByNonResidue(x *E4) *E4 {
+	z.B1, z.B0 = x.B0, x.B1
+	z.B0.MulByNonResidue(&z.B0)
+	return z
+}
+
+// Mul sets z=x*y in E4 and returns z
+func (z *E4) Mul(x, y *E4) *E4 {
+	var a, b, c E2
+	a.Add(&x.B0, &x.B1)
+	b.Add(&y.B0, &y.B1)
+	a.Mul(&a, &b)
+	b.Mul(&x.B0, &y.B0)
+	c.Mul(&x.B1, &y.B1)
+	z.B1.Sub(&a, &b).Sub(&z.B1, &c)
+	z.B0.MulByNonResidue(&c).Add(&z.B0, &b)
+	return z
+}
+
+// Square sets z=x*x in E4 and returns z
+func (z *E4) Square(x *E4) *E4 {
+
+	//Algorithm 22 from https://eprint.iacr.org/2010/354.pdf
+	var c0, c2, c3 E2
+	c0.Sub(&x.B0, &x.B1)
+	c3.MulByNonResidue(&x.B1).Sub(&x.B0, &c3)
+	c2.Mul(&x.B0, &x.B1)
+	c0.Mul(&c0, &c3).Add(&c0, &c2)
+	z.B1.Double(&c2)
+	c2.MulByNonResidue(&c2)
+	z.B0.Add(&c0, &c2)
+
+	return z
+}
+
+// Inverse sets z to the inverse of x in E4 and returns z
+//
+// if x == 0, sets and returns z = x
+func (z *E4) Inverse(x *E4) *E4 {
+	// Algorithm 23 from https://eprint.iacr.org/2010/354.pdf
+
+	var t0, t1, tmp E2
+	t0.Square(&x.B0)
+	t1.Square(&x.B1)
+	tmp.MulByNonResidue(&t1)
+	t0.Sub(&t0, &tmp)
+	t1.Inverse(&t0)
+	z.B0.Mul(&x.B0, &t1)
+	z.B1.Mul(&x.B1, &t1).Neg(&z.B1)
+
+	return z
+}
+
+// Exp sets z=xᵏ (mod q⁴) and returns it
+func (z *E4) Exp(x E4, k *big.Int) *E4 {
+	if k.IsUint64() && k.Uint64() == 0 {
+		return z.SetOne()
+	}
+
+	e := k
+	if k.Sign() == -1 {
+		// negative k, we invert
+		// if k < 0: xᵏ (mod q⁴) == (x⁻¹)ᵏ (mod q⁴)
+		x.Inverse(&x)
+
+		// we negate k in a temp big.Int since
+		// Int.Bit(_) of k and -k is different
+		e = bigIntPool.Get().(*big.Int)
+		defer bigIntPool.Put(e)
+		e.Neg(k)
+	}
+
+	z.SetOne()
+	b := e.Bytes()
+	for i := 0; i < len(b); i++ {
+		w := b[i]
+		for j := 0; j < 8; j++ {
+			z.Square(z)
+			if (w & (0b10000000 >> j)) != 0 {
+				z.Mul(z, &x)
+			}
+		}
+	}
+
+	return z
+}
+
+// Conjugate sets z to x conjugated and returns z
+func (z *E4) Conjugate(x *E4) *E4 {
+	z.B0 = x.B0
+	z.B1.Neg(&x.B1)
+	return z
+}
+
+func (z *E4) Halve() {
+
+	z.B0.A0.Halve()
+	z.B0.A1.Halve()
+	z.B1.A0.Halve()
+	z.B1.A1.Halve()
+}
+
+// norm sets x to the norm of z
+func (z *E4) norm(x *E2) {
+	var tmp E2
+	tmp.Square(&z.B1).MulByNonResidue(&tmp)
+	x.Square(&z.B0).Sub(x, &tmp)
+}
+
+// Legendre returns the Legendre symbol of z
+func (z *E4) Legendre() int {
+	var n E2
+	z.norm(&n)
+	return n.Legendre()
+}
+
+// Sqrt sets z to the square root of and returns z
+// The function does not test whether the square root
+// exists or not, it's up to the caller to call
+// Legendre beforehand.
+// cf https://eprint.iacr.org/2012/685.pdf (algo 10)
+func (z *E4) Sqrt(x *E4) *E4 {
+
+	// precomputation
+	var b, c, d, e, f, x0, _g E4
+	var _b, o E2
+
+	// c must be a non square (p = 1 mod 4)
+	c.B1.SetOne()
+
+	q := fr.Modulus()
+	var exp, one big.Int
+	one.SetUint64(1)
+	exp.Mul(q, q).Sub(&exp, &one).Rsh(&exp, 1)
+	d.Exp(c, &exp)
+	e.Mul(&d, &c).Inverse(&e)
+	f.Mul(&d, &c).Square(&f)
+
+	// computation
+	exp.Rsh(&exp, 1)
+	b.Exp(*x, &exp)
+	b.norm(&_b)
+	o.SetOne()
+	if _b.Equal(&o) {
+		x0.Square(&b).Mul(&x0, x)
+		_b.Set(&x0.B0).Sqrt(&_b)
+		_g.B0.Set(&_b)
+		z.Conjugate(&b).Mul(z, &_g)
+		return z
+	}
+	x0.Square(&b).Mul(&x0, x).Mul(&x0, &f)
+	_b.Set(&x0.B0).Sqrt(&_b)
+	_g.B0.Set(&_b)
+	z.Conjugate(&b).Mul(z, &_g).Mul(z, &e)
+
+	return z
+}
+
+// BatchInvertE4 returns a new slice with every element in a inverted.
+// It uses Montgomery batch inversion trick.
+//
+// if a[i] == 0, returns result[i] = a[i]
+func BatchInvertE4(a []E4) []E4 {
+	res := make([]E4, len(a))
+	if len(a) == 0 {
+		return res
+	}
+
+	zeroes := make([]bool, len(a))
+	var accumulator E4
+	accumulator.SetOne()
+
+	for i := 0; i < len(a); i++ {
+		if a[i].IsZero() {
+			zeroes[i] = true
+			continue
+		}
+		res[i].Set(&accumulator)
+		accumulator.Mul(&accumulator, &a[i])
+	}
+
+	accumulator.Inverse(&accumulator)
+
+	for i := len(a) - 1; i >= 0; i-- {
+		if zeroes[i] {
+			continue
+		}
+		res[i].Mul(&res[i], &accumulator)
+		accumulator.Mul(&accumulator, &a[i])
+	}
+
+	return res
+}
+
+// Div divides an element in E4 by an element in E4
+func (z *E4) Div(x *E4, y *E4) *E4 {
+	var r E4
+	r.Inverse(y).Mul(x, &r)
+	return z.Set(&r)
+}
diff --git a/field/babybear/extensions/e4_test.go b/field/babybear/extensions/e4_test.go
new file mode 100644
index 0000000000..69d43787f2
--- /dev/null
+++ b/field/babybear/extensions/e4_test.go
@@ -0,0 +1,300 @@
+package extensions
+
+import (
+	"testing"
+
+	"github.com/leanovate/gopter"
+	"github.com/leanovate/gopter/prop"
+)
+
+// ------------------------------------------------------------
+// tests
+
+func TestE4ReceiverIsOperand(t *testing.T) {
+
+	parameters := gopter.DefaultTestParameters()
+	parameters.MinSuccessfulTests = 100
+
+	properties := gopter.NewProperties(parameters)
+
+	genA := GenE4()
+	genB := GenE4()
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (addition) should output the same result", prop.ForAll(
+		func(a, b *E4) bool {
+			var c, d E4
+			d.Set(a)
+			c.Add(a, b)
+			a.Add(a, b)
+			b.Add(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (sub) should output the same result", prop.ForAll(
+		func(a, b *E4) bool {
+			var c, d E4
+			d.Set(a)
+			c.Sub(a, b)
+			a.Sub(a, b)
+			b.Sub(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (mul) should output the same result", prop.ForAll(
+		func(a, b *E4) bool {
+			var c, d E4
+			d.Set(a)
+			c.Mul(a, b)
+			a.Mul(a, b)
+			b.Mul(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (square) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Square(a)
+			a.Square(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (double) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Double(a)
+			a.Double(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (mul by non residue) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.MulByNonResidue(a)
+			a.MulByNonResidue(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (Inverse) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Inverse(a)
+			a.Inverse(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (Conjugate) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Conjugate(a)
+			a.Conjugate(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Having the receiver as operand (Sqrt) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b, c, d, s E4
+
+			s.Square(a)
+			a.Set(&s)
+			b.Set(&s)
+
+			a.Sqrt(a)
+			b.Sqrt(&b)
+
+			c.Square(a)
+			d.Square(&b)
+			return c.Equal(&d)
+		},
+		genA,
+	))
+
+	properties.TestingRun(t, gopter.ConsoleReporter(false))
+}
+
+func TestE4Ops(t *testing.T) {
+
+	parameters := gopter.DefaultTestParameters()
+	parameters.MinSuccessfulTests = 100
+
+	properties := gopter.NewProperties(parameters)
+
+	genA := GenE4()
+	genB := GenE4()
+
+	properties.Property("[BABYBEAR] sub & add should leave an element invariant", prop.ForAll(
+		func(a, b *E4) bool {
+			var c E4
+			c.Set(a)
+			c.Add(&c, b).Sub(&c, b)
+			return c.Equal(a)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[BABYBEAR] mul & inverse should leave an element invariant", prop.ForAll(
+		func(a, b *E4) bool {
+			var c, d E4
+			d.Inverse(b)
+			c.Set(a)
+			c.Mul(&c, b).Mul(&c, &d)
+			return c.Equal(a)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[BABYBEAR] BatchInvertE4 should output the same result as Inverse", prop.ForAll(
+		func(a, b, c *E4) bool {
+
+			batch := BatchInvertE4([]E4{*a, *b, *c})
+			a.Inverse(a)
+			b.Inverse(b)
+			c.Inverse(c)
+			return a.Equal(&batch[0]) && b.Equal(&batch[1]) && c.Equal(&batch[2])
+		},
+		genA,
+		genA,
+		genB,
+	))
+
+	properties.Property("[BABYBEAR] inverse twice should leave an element invariant", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Inverse(a).Inverse(&b)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] square and mul should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b, c E4
+			b.Mul(a, a)
+			c.Square(a)
+			return b.Equal(&c)
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] Legendre on square should output 1", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Square(a)
+			c := b.Legendre()
+			return c == 1
+		},
+		genA,
+	))
+
+	properties.Property("[BABYBEAR] square(sqrt) should leave an element invariant", prop.ForAll(
+		func(a *E4) bool {
+			var b, c, d, e E4
+			b.Square(a)
+			c.Sqrt(&b)
+			d.Square(&c)
+			e.Neg(a)
+			return (c.Equal(a) || c.Equal(&e)) && d.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.TestingRun(t, gopter.ConsoleReporter(false))
+}
+
+// ------------------------------------------------------------
+// benches
+
+func BenchmarkE4Add(b *testing.B) {
+	var a, c E4
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Add(&a, &c)
+	}
+}
+
+func BenchmarkE4Sub(b *testing.B) {
+	var a, c E4
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Sub(&a, &c)
+	}
+}
+
+func BenchmarkE4Mul(b *testing.B) {
+	var a, c E4
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Mul(&a, &c)
+	}
+}
+
+func BenchmarkE4Square(b *testing.B) {
+	var a E4
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Square(&a)
+	}
+}
+
+func BenchmarkE4Sqrt(b *testing.B) {
+	var a E4
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Sqrt(&a)
+	}
+}
+
+func BenchmarkE4Inverse(b *testing.B) {
+	var a E4
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Inverse(&a)
+	}
+}
+
+func BenchmarkE4MulNonRes(b *testing.B) {
+	var a E4
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.MulByNonResidue(&a)
+	}
+}
+
+func BenchmarkE4Conjugate(b *testing.B) {
+	var a E4
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Conjugate(&a)
+	}
+}
diff --git a/field/babybear/extensions/utils.go b/field/babybear/extensions/utils.go
new file mode 100644
index 0000000000..1cb2425e1f
--- /dev/null
+++ b/field/babybear/extensions/utils.go
@@ -0,0 +1,48 @@
+package extensions
+
+import (
+	"math/big"
+	"sync"
+
+	fr "github.com/consensys/gnark-crypto/field/babybear"
+	"github.com/leanovate/gopter"
+)
+
+var bigIntPool = sync.Pool{
+	New: func() interface{} {
+		return new(big.Int)
+	},
+}
+
+// Fr generates an Fr element
+func GenFr() gopter.Gen {
+	return func(genParams *gopter.GenParameters) *gopter.GenResult {
+		var elmt fr.Element
+
+		if _, err := elmt.SetRandom(); err != nil {
+			panic(err)
+		}
+		genResult := gopter.NewGenResult(elmt, gopter.NoShrinker)
+		return genResult
+	}
+}
+
+// E2 generates an E2 elmt
+func GenE2() gopter.Gen {
+	return gopter.CombineGens(
+		GenFr(),
+		GenFr(),
+	).Map(func(values []interface{}) *E2 {
+		return &E2{A0: values[0].(fr.Element), A1: values[1].(fr.Element)}
+	})
+}
+
+// E4 generates an E4 elmt
+func GenE4() gopter.Gen {
+	return gopter.CombineGens(
+		GenE2(),
+		GenE2(),
+	).Map(func(values []interface{}) *E4 {
+		return &E4{B0: *values[0].(*E2), B1: *values[1].(*E2)}
+	})
+}
diff --git a/field/generator/config/field_config.go b/field/generator/config/field_config.go
index 5cbc09dac8..a8485cb8e3 100644
--- a/field/generator/config/field_config.go
+++ b/field/generator/config/field_config.go
@@ -42,7 +42,7 @@ type Field struct {
 	QMinusOneHalvedP          []uint64 // ((q-1) / 2 ) + 1
 	Mu                        uint64   // mu = 2^288 / q for 4.5 word barrett reduction
 	RSquare                   []uint64
-	One, Thirteen             []uint64
+	One, Eleven, Thirteen     []uint64
 	LegendreExponent          string // big.Int to base16 string
 	NoCarry                   bool
 	NoCarrySquare             bool // used if NoCarry is set, but some op may overflow in square optimization
@@ -182,6 +182,13 @@ func NewFieldConfig(packageName, elementName, modulus string, useAddChain bool)
 	one.Lsh(&one, uint(F.NbWords)*radix).Mod(&one, &bModulus)
 	F.One = toUint64Slice(&one, F.NbWords)
 
+	{
+		var n big.Int
+		n.SetUint64(11)
+		n.Lsh(&n, uint(F.NbWords)*radix).Mod(&n, &bModulus)
+		F.Eleven = toUint64Slice(&n, F.NbWords)
+	}
+
 	{
 		var n big.Int
 		n.SetUint64(13)
diff --git a/field/koalabear/extensions/doc.go b/field/koalabear/extensions/doc.go
new file mode 100644
index 0000000000..f3a8775be3
--- /dev/null
+++ b/field/koalabear/extensions/doc.go
@@ -0,0 +1,6 @@
+// Package extensions implements the fields arithmetic of the 𝔽r² and 𝔽r⁴
+// extensions of the Koala-bear field.
+//
+//	𝔽r²[u] = 𝔽r/u²-4
+//	𝔽r⁴[v] = 𝔽r²/v²-u
+package extensions
diff --git a/field/koalabear/extensions/e2.go b/field/koalabear/extensions/e2.go
new file mode 100644
index 0000000000..5c6df7b218
--- /dev/null
+++ b/field/koalabear/extensions/e2.go
@@ -0,0 +1,278 @@
+package extensions
+
+import (
+	"math/big"
+
+	fr "github.com/consensys/gnark-crypto/field/koalabear"
+)
+
+// E2 is a degree two finite field extension of fr.Element
+type E2 struct {
+	A0, A1 fr.Element
+}
+
+// Equal returns true if z equals x, false otherwise
+func (z *E2) Equal(x *E2) bool {
+	return z.A0.Equal(&x.A0) && z.A1.Equal(&x.A1)
+}
+
+// Cmp compares (lexicographic order) z and x and returns:
+//
+//	-1 if z <  x
+//	 0 if z == x
+//	+1 if z >  x
+func (z *E2) Cmp(x *E2) int {
+	if a1 := z.A1.Cmp(&x.A1); a1 != 0 {
+		return a1
+	}
+	return z.A0.Cmp(&x.A0)
+}
+
+// LexicographicallyLargest returns true if this element is strictly lexicographically
+// larger than its negation, false otherwise
+func (z *E2) LexicographicallyLargest() bool {
+	// adapted from github.com/zkcrypto/bls12_381
+	if z.A1.IsZero() {
+		return z.A0.LexicographicallyLargest()
+	}
+	return z.A1.LexicographicallyLargest()
+}
+
+// SetString sets a E2 element from strings
+func (z *E2) SetString(s1, s2 string) *E2 {
+	z.A0.SetString(s1)
+	z.A1.SetString(s2)
+	return z
+}
+
+// SetZero sets an E2 elmt to zero
+func (z *E2) SetZero() *E2 {
+	z.A0.SetZero()
+	z.A1.SetZero()
+	return z
+}
+
+// Set sets an E2 from x
+func (z *E2) Set(x *E2) *E2 {
+	z.A0 = x.A0
+	z.A1 = x.A1
+	return z
+}
+
+// SetOne sets z to 1 in Montgomery form and returns z
+func (z *E2) SetOne() *E2 {
+	z.A0.SetOne()
+	z.A1.SetZero()
+	return z
+}
+
+// SetRandom sets a0 and a1 to random values
+func (z *E2) SetRandom() (*E2, error) {
+	if _, err := z.A0.SetRandom(); err != nil {
+		return nil, err
+	}
+	if _, err := z.A1.SetRandom(); err != nil {
+		return nil, err
+	}
+	return z, nil
+}
+
+// IsZero returns true if z is zero, false otherwise
+func (z *E2) IsZero() bool {
+	return z.A0.IsZero() && z.A1.IsZero()
+}
+
+// IsOne returns true if z is one, false otherwise
+func (z *E2) IsOne() bool {
+	return z.A0.IsOne() && z.A1.IsZero()
+}
+
+// Add adds two elements of E2
+func (z *E2) Add(x, y *E2) *E2 {
+	z.A0.Add(&x.A0, &y.A0)
+	z.A1.Add(&x.A1, &y.A1)
+	return z
+}
+
+// Sub subtracts two elements of E2
+func (z *E2) Sub(x, y *E2) *E2 {
+	z.A0.Sub(&x.A0, &y.A0)
+	z.A1.Sub(&x.A1, &y.A1)
+	return z
+}
+
+// Double doubles an E2 element
+func (z *E2) Double(x *E2) *E2 {
+	z.A0.Double(&x.A0)
+	z.A1.Double(&x.A1)
+	return z
+}
+
+// Neg negates an E2 element
+func (z *E2) Neg(x *E2) *E2 {
+	z.A0.Neg(&x.A0)
+	z.A1.Neg(&x.A1)
+	return z
+}
+
+// String implements Stringer interface for fancy printing
+func (z *E2) String() string {
+	return z.A0.String() + "+" + z.A1.String() + "*u"
+}
+
+// MulByElement multiplies an element in E2 by an element in fr
+func (z *E2) MulByElement(x *E2, y *fr.Element) *E2 {
+	var yCopy fr.Element
+	yCopy.Set(y)
+	z.A0.Mul(&x.A0, &yCopy)
+	z.A1.Mul(&x.A1, &yCopy)
+	return z
+}
+
+// Conjugate conjugates an element in E2
+func (z *E2) Conjugate(x *E2) *E2 {
+	z.A0 = x.A0
+	z.A1.Neg(&x.A1)
+	return z
+}
+
+// Halve sets z to z / 2
+func (z *E2) Halve() {
+	z.A0.Halve()
+	z.A1.Halve()
+}
+
+// Legendre returns the Legendre symbol of z
+func (z *E2) Legendre() int {
+	var n fr.Element
+	z.norm(&n)
+	return n.Legendre()
+}
+
+// Exp sets z=xᵏ (mod q²) and returns it
+func (z *E2) Exp(x E2, k *big.Int) *E2 {
+	if k.IsUint64() && k.Uint64() == 0 {
+		return z.SetOne()
+	}
+
+	e := k
+	if k.Sign() == -1 {
+		// negative k, we invert
+		// if k < 0: xᵏ (mod q²) == (x⁻¹)ᵏ (mod q²)
+		x.Inverse(&x)
+
+		// we negate k in a temp big.Int since
+		// Int.Bit(_) of k and -k is different
+		e = bigIntPool.Get().(*big.Int)
+		defer bigIntPool.Put(e)
+		e.Neg(k)
+	}
+
+	z.SetOne()
+	b := e.Bytes()
+	for i := 0; i < len(b); i++ {
+		w := b[i]
+		for j := 0; j < 8; j++ {
+			z.Square(z)
+			if (w & (0b10000000 >> j)) != 0 {
+				z.Mul(z, &x)
+			}
+		}
+	}
+
+	return z
+}
+
+// Sqrt sets z to the square root of and returns z
+// The function does not test whether the square root
+// exists or not, it's up to the caller to call
+// Legendre beforehand.
+// cf https://eprint.iacr.org/2012/685.pdf (algo 10)
+func (z *E2) Sqrt(x *E2) *E2 {
+
+	// precomputation
+	var b, c, d, e, f, x0 E2
+	var _b, o fr.Element
+
+	// c must be a non square (p = 1 mod 4)
+	c.A1.SetOne()
+
+	q := fr.Modulus()
+	var exp, one big.Int
+	one.SetUint64(1)
+	exp.Set(q).Sub(&exp, &one).Rsh(&exp, 1)
+	d.Exp(c, &exp)
+	e.Mul(&d, &c).Inverse(&e)
+	f.Mul(&d, &c).Square(&f)
+
+	// computation
+	exp.Rsh(&exp, 1)
+	b.Exp(*x, &exp)
+	b.norm(&_b)
+	o.SetOne()
+	if _b.Equal(&o) {
+		x0.Square(&b).Mul(&x0, x)
+		_b.Set(&x0.A0).Sqrt(&_b)
+		z.Conjugate(&b).MulByElement(z, &_b)
+		return z
+	}
+	x0.Square(&b).Mul(&x0, x).Mul(&x0, &f)
+	_b.Set(&x0.A0).Sqrt(&_b)
+	z.Conjugate(&b).MulByElement(z, &_b).Mul(z, &e)
+
+	return z
+}
+
+// BatchInvertE2 returns a new slice with every element in a inverted.
+// It uses Montgomery batch inversion trick.
+//
+// if a[i] == 0, returns result[i] = a[i]
+func BatchInvertE2(a []E2) []E2 {
+	res := make([]E2, len(a))
+	if len(a) == 0 {
+		return res
+	}
+
+	zeroes := make([]bool, len(a))
+	var accumulator E2
+	accumulator.SetOne()
+
+	for i := 0; i < len(a); i++ {
+		if a[i].IsZero() {
+			zeroes[i] = true
+			continue
+		}
+		res[i].Set(&accumulator)
+		accumulator.Mul(&accumulator, &a[i])
+	}
+
+	accumulator.Inverse(&accumulator)
+
+	for i := len(a) - 1; i >= 0; i-- {
+		if zeroes[i] {
+			continue
+		}
+		res[i].Mul(&res[i], &accumulator)
+		accumulator.Mul(&accumulator, &a[i])
+	}
+
+	return res
+}
+
+// Select is conditional move.
+// If cond = 0, it sets z to caseZ and returns it. otherwise caseNz.
+func (z *E2) Select(cond int, caseZ *E2, caseNz *E2) *E2 {
+	//Might be able to save a nanosecond or two by an aggregate implementation
+
+	z.A0.Select(cond, &caseZ.A0, &caseNz.A0)
+	z.A1.Select(cond, &caseZ.A1, &caseNz.A1)
+
+	return z
+}
+
+// Div divides an element in E2 by an element in E2
+func (z *E2) Div(x *E2, y *E2) *E2 {
+	var r E2
+	r.Inverse(y).Mul(x, &r)
+	return z.Set(&r)
+}
diff --git a/field/koalabear/extensions/e2_koalabear.go b/field/koalabear/extensions/e2_koalabear.go
new file mode 100644
index 0000000000..68d04c7a08
--- /dev/null
+++ b/field/koalabear/extensions/e2_koalabear.go
@@ -0,0 +1,79 @@
+package extensions
+
+import (
+	fr "github.com/consensys/gnark-crypto/field/koalabear"
+)
+
+// Mul sets z to the E2-product of x,y, returns z
+func (z *E2) Mul(x, y *E2) *E2 {
+	var a, b, c fr.Element
+	a.Add(&x.A0, &x.A1)
+	b.Add(&y.A0, &y.A1)
+	a.Mul(&a, &b)
+	b.Mul(&x.A0, &y.A0)
+	c.Mul(&x.A1, &y.A1)
+	z.A1.Sub(&a, &b).Sub(&z.A1, &c)
+	fr.MulBy3(&c)
+	z.A0.Add(&b, &c)
+	return z
+}
+
+// Square sets z to the E2-product of x,x returns z
+func (z *E2) Square(x *E2) *E2 {
+	var a, b, c fr.Element
+	a.Mul(&x.A0, &x.A1).Double(&a)
+	c.Square(&x.A0)
+	b.Square(&x.A1)
+	fr.MulBy3(&b)
+	z.A0.Add(&c, &b)
+	z.A1 = a
+	return z
+}
+
+// MulByNonResidue multiplies a E2 by (0,1)
+func (z *E2) MulByNonResidue(x *E2) *E2 {
+	a := x.A0
+	b := x.A1 // fetching x.A1 in the function below is slower
+	fr.MulBy3(&b)
+	z.A0 = b
+	z.A1 = a
+	return z
+}
+
+// MulByNonResidueInv multiplies a E2 by (0,1)^{-1}
+func (z *E2) MulByNonResidueInv(x *E2) *E2 {
+	a := x.A1
+	// 1/3 mod r
+	var threeInv fr.Element
+	threeInv.SetUint64(710235478)
+	z.A1.Mul(&x.A0, &threeInv)
+	z.A0 = a
+	return z
+}
+
+// Inverse sets z to the E2-inverse of x, returns z
+func (z *E2) Inverse(x *E2) *E2 {
+	// Algorithm 8 from https://eprint.iacr.org/2010/354.pdf
+	var t0, t1, tmp fr.Element
+	a := &x.A0 // creating the buffers a, b is faster than querying &x.A0, &x.A1 in the functions call below
+	b := &x.A1
+	t0.Square(a)
+	t1.Square(b)
+	tmp.Set(&t1)
+	fr.MulBy3(&tmp)
+	t0.Sub(&t0, &tmp)
+	t1.Inverse(&t0)
+	z.A0.Mul(a, &t1)
+	z.A1.Mul(b, &t1).Neg(&z.A1)
+
+	return z
+}
+
+// norm sets x to the norm of z
+func (z *E2) norm(x *fr.Element) {
+	var tmp fr.Element
+	x.Square(&z.A1)
+	tmp.Set(x)
+	fr.MulBy3(&tmp)
+	x.Square(&z.A0).Sub(x, &tmp)
+}
diff --git a/field/koalabear/extensions/e2_test.go b/field/koalabear/extensions/e2_test.go
new file mode 100644
index 0000000000..998b5e6a51
--- /dev/null
+++ b/field/koalabear/extensions/e2_test.go
@@ -0,0 +1,510 @@
+package extensions
+
+import (
+	"crypto/rand"
+	"testing"
+
+	fr "github.com/consensys/gnark-crypto/field/koalabear"
+	"github.com/leanovate/gopter"
+	"github.com/leanovate/gopter/prop"
+)
+
+// ------------------------------------------------------------
+// tests
+
+const (
+	nbFuzzShort = 10
+	nbFuzz      = 50
+)
+
+func TestE2ReceiverIsOperand(t *testing.T) {
+
+	t.Parallel()
+	parameters := gopter.DefaultTestParameters()
+	if testing.Short() {
+		parameters.MinSuccessfulTests = nbFuzzShort
+	} else {
+		parameters.MinSuccessfulTests = nbFuzz
+	}
+
+	properties := gopter.NewProperties(parameters)
+
+	genA := GenE2()
+	genB := GenE2()
+	genfr := GenFr()
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (addition) should output the same result", prop.ForAll(
+		func(a, b *E2) bool {
+			var c, d E2
+			d.Set(a)
+			c.Add(a, b)
+			a.Add(a, b)
+			b.Add(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (sub) should output the same result", prop.ForAll(
+		func(a, b *E2) bool {
+			var c, d E2
+			d.Set(a)
+			c.Sub(a, b)
+			a.Sub(a, b)
+			b.Sub(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (mul) should output the same result", prop.ForAll(
+		func(a, b *E2) bool {
+			var c, d E2
+			d.Set(a)
+			c.Mul(a, b)
+			a.Mul(a, b)
+			b.Mul(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (square) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Square(a)
+			a.Square(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (neg) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Neg(a)
+			a.Neg(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (double) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Double(a)
+			a.Double(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (mul by non residue) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.MulByNonResidue(a)
+			a.MulByNonResidue(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (mul by non residue inverse) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.MulByNonResidueInv(a)
+			a.MulByNonResidueInv(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (Inverse) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Inverse(a)
+			a.Inverse(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (Conjugate) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Conjugate(a)
+			a.Conjugate(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (mul by element) should output the same result", prop.ForAll(
+		func(a *E2, b fr.Element) bool {
+			var c E2
+			c.MulByElement(a, &b)
+			a.MulByElement(a, &b)
+			return a.Equal(&c)
+		},
+		genA,
+		genfr,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (Sqrt) should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b, c, d, s E2
+
+			s.Square(a)
+			a.Set(&s)
+			b.Set(&s)
+
+			a.Sqrt(a)
+			b.Sqrt(&b)
+
+			c.Square(a)
+			d.Square(&b)
+			return c.Equal(&d)
+		},
+		genA,
+	))
+
+	properties.TestingRun(t, gopter.ConsoleReporter(false))
+
+}
+
+func TestE2MulMaxed(t *testing.T) {
+	// let's pick a and b, with maxed A0 and A1
+	var a, b E2
+	frMaxValue := fr.Element{
+		2130706433,
+	}
+	frMaxValue[0]--
+
+	a.A0 = frMaxValue
+	a.A1 = frMaxValue
+	b.A0 = frMaxValue
+	b.A1 = frMaxValue
+
+	var c, d E2
+	d.Inverse(&b)
+	c.Set(&a)
+	c.Mul(&c, &b).Mul(&c, &d)
+	if !c.Equal(&a) {
+		t.Fatal("mul with max fr failed")
+	}
+}
+
+func TestE2Ops(t *testing.T) {
+
+	t.Parallel()
+	parameters := gopter.DefaultTestParameters()
+	if testing.Short() {
+		parameters.MinSuccessfulTests = nbFuzzShort
+	} else {
+		parameters.MinSuccessfulTests = nbFuzz
+	}
+
+	properties := gopter.NewProperties(parameters)
+
+	genA := GenE2()
+	genB := GenE2()
+	genfr := GenFr()
+
+	properties.Property("[KOALABEAR] sub & add should leave an element invariant", prop.ForAll(
+		func(a, b *E2) bool {
+			var c E2
+			c.Set(a)
+			c.Add(&c, b).Sub(&c, b)
+			return c.Equal(a)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[KOALABEAR] mul & inverse should leave an element invariant", prop.ForAll(
+		func(a, b *E2) bool {
+			var c, d E2
+			d.Inverse(b)
+			c.Set(a)
+			c.Mul(&c, b).Mul(&c, &d)
+			return c.Equal(a)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[KOALABEAR] BatchInvertE2 should output the same result as Inverse", prop.ForAll(
+		func(a, b, c *E2) bool {
+
+			batch := BatchInvertE2([]E2{*a, *b, *c})
+			a.Inverse(a)
+			b.Inverse(b)
+			c.Inverse(c)
+			return a.Equal(&batch[0]) && b.Equal(&batch[1]) && c.Equal(&batch[2])
+		},
+		genA,
+		genA,
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] inverse twice should leave an element invariant", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Inverse(a).Inverse(&b)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] neg twice should leave an element invariant", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Neg(a).Neg(&b)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] square and mul should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b, c E2
+			b.Mul(a, a)
+			c.Square(a)
+			return b.Equal(&c)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] MulByElement MulByElement inverse should leave an element invariant", prop.ForAll(
+		func(a *E2, b fr.Element) bool {
+			var c E2
+			var d fr.Element
+			d.Inverse(&b)
+			c.MulByElement(a, &b).MulByElement(&c, &d)
+			return c.Equal(a)
+		},
+		genA,
+		genfr,
+	))
+
+	properties.Property("[KOALABEAR] Double and mul by 2 should output the same result", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			var c fr.Element
+			c.SetUint64(2)
+			b.Double(a)
+			a.MulByElement(a, &c)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Mulbynonres mulbynonresinv should leave the element invariant", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.MulByNonResidue(a).MulByNonResidueInv(&b)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] a + pi(a), a-pi(a) should be real", prop.ForAll(
+		func(a *E2) bool {
+			var b, c, d E2
+			var e, f fr.Element
+			b.Conjugate(a)
+			c.Add(a, &b)
+			d.Sub(a, &b)
+			e.Double(&a.A0)
+			f.Double(&a.A1)
+			return c.A1.IsZero() && d.A0.IsZero() && e.Equal(&c.A0) && f.Equal(&d.A1)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Legendre on square should output 1", prop.ForAll(
+		func(a *E2) bool {
+			var b E2
+			b.Square(a)
+			c := b.Legendre()
+			return c == 1
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] square(sqrt) should leave an element invariant", prop.ForAll(
+		func(a *E2) bool {
+			var b, c, d, e E2
+			b.Square(a)
+			c.Sqrt(&b)
+			d.Square(&c)
+			e.Neg(a)
+			return (c.Equal(a) || c.Equal(&e)) && d.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] neg(E2) == neg(E2.A0, E2.A1)", prop.ForAll(
+		func(a *E2) bool {
+			var b, c E2
+			b.Neg(a)
+			c.A0.Neg(&a.A0)
+			c.A1.Neg(&a.A1)
+			return c.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Cmp and LexicographicallyLargest should be consistent", prop.ForAll(
+		func(a *E2) bool {
+			var negA E2
+			negA.Neg(a)
+			cmpResult := a.Cmp(&negA)
+			lResult := a.LexicographicallyLargest()
+			if lResult && cmpResult == 1 {
+				return true
+			}
+			if !lResult && cmpResult != 1 {
+				return true
+			}
+			return false
+		},
+		genA,
+	))
+
+	properties.TestingRun(t, gopter.ConsoleReporter(false))
+
+}
+
+// ------------------------------------------------------------
+// benches
+
+func BenchmarkE2Add(b *testing.B) {
+	var a, c E2
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Add(&a, &c)
+	}
+}
+
+func BenchmarkE2Sub(b *testing.B) {
+	var a, c E2
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Sub(&a, &c)
+	}
+}
+
+func BenchmarkE2Mul(b *testing.B) {
+	var a, c E2
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Mul(&a, &c)
+	}
+}
+
+func BenchmarkE2MulByElement(b *testing.B) {
+	var a E2
+	var c fr.Element
+	_, _ = c.SetRandom()
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.MulByElement(&a, &c)
+	}
+}
+
+func BenchmarkE2Square(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Square(&a)
+	}
+}
+
+func BenchmarkE2Sqrt(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Sqrt(&a)
+	}
+}
+
+func BenchmarkE2Exp(b *testing.B) {
+	var x E2
+	_, _ = x.SetRandom()
+	b1, _ := rand.Int(rand.Reader, fr.Modulus())
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		x.Exp(x, b1)
+	}
+}
+
+func BenchmarkE2Inverse(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Inverse(&a)
+	}
+}
+
+func BenchmarkE2MulNonRes(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.MulByNonResidue(&a)
+	}
+}
+
+func BenchmarkE2MulNonResInv(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.MulByNonResidueInv(&a)
+	}
+}
+
+func BenchmarkE2Conjugate(b *testing.B) {
+	var a E2
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Conjugate(&a)
+	}
+}
+
+func TestE2Div(t *testing.T) {
+
+	parameters := gopter.DefaultTestParameters()
+	properties := gopter.NewProperties(parameters)
+
+	genA := GenE2()
+	genB := GenE2()
+
+	properties.Property("[KOALABEAR] dividing then multiplying by the same element does nothing", prop.ForAll(
+		func(a, b *E2) bool {
+			var c E2
+			c.Div(a, b)
+			c.Mul(&c, b)
+			return c.Equal(a)
+		},
+		genA,
+		genB,
+	))
+
+	properties.TestingRun(t, gopter.ConsoleReporter(false))
+}
diff --git a/field/koalabear/extensions/e4.go b/field/koalabear/extensions/e4.go
new file mode 100644
index 0000000000..c6fd21bd09
--- /dev/null
+++ b/field/koalabear/extensions/e4.go
@@ -0,0 +1,341 @@
+package extensions
+
+import (
+	"math/big"
+
+	fr "github.com/consensys/gnark-crypto/field/koalabear"
+)
+
+// E4 is a degree two finite field extension of fr2
+type E4 struct {
+	B0, B1 E2
+}
+
+// Equal returns true if z equals x, false otherwise
+func (z *E4) Equal(x *E4) bool {
+	return z.B0.Equal(&x.B0) && z.B1.Equal(&x.B1)
+}
+
+// Cmp compares (lexicographic order) z and x and returns:
+//
+//	-1 if z <  x
+//	 0 if z == x
+//	+1 if z >  x
+func (z *E4) Cmp(x *E4) int {
+	if a1 := z.B1.Cmp(&x.B1); a1 != 0 {
+		return a1
+	}
+	return z.B0.Cmp(&x.B0)
+}
+
+// LexicographicallyLargest returns true if this element is strictly lexicographically
+// larger than its negation, false otherwise
+func (z *E4) LexicographicallyLargest() bool {
+	// adapted from github.com/zkcrypto/bls12_381
+	if z.B1.IsZero() {
+		return z.B0.LexicographicallyLargest()
+	}
+	return z.B1.LexicographicallyLargest()
+}
+
+// String puts E4 in string form
+func (z *E4) String() string {
+	return (z.B0.String() + "+(" + z.B1.String() + ")*v")
+}
+
+// SetString sets a E4 from string
+func (z *E4) SetString(s0, s1, s2, s3 string) *E4 {
+	z.B0.SetString(s0, s1)
+	z.B1.SetString(s2, s3)
+	return z
+}
+
+// Set copies x into z and returns z
+func (z *E4) Set(x *E4) *E4 {
+	z.B0 = x.B0
+	z.B1 = x.B1
+	return z
+}
+
+// SetZero sets an E4 elmt to zero
+func (z *E4) SetZero() *E4 {
+	z.B0.SetZero()
+	z.B1.SetZero()
+	return z
+}
+
+// SetOne sets z to 1 in Montgomery form and returns z
+func (z *E4) SetOne() *E4 {
+	*z = E4{}
+	z.B0.A0.SetOne()
+	return z
+}
+
+// MulByElement multiplies an element in E4 by an element in fr
+func (z *E4) MulByElement(x *E4, y *fr.Element) *E4 {
+	var yCopy fr.Element
+	yCopy.Set(y)
+	z.B0.MulByElement(&x.B0, &yCopy)
+	z.B1.MulByElement(&x.B1, &yCopy)
+	return z
+}
+
+// MulByE2 multiplies an element in E4 by an element in E2
+func (z *E4) MulByE2(x *E4, y *E2) *E4 {
+	var yCopy E2
+	yCopy.Set(y)
+	z.B0.Mul(&x.B0, &yCopy)
+	z.B1.Mul(&x.B1, &yCopy)
+	return z
+}
+
+// Add sets z=x+y in E4 and returns z
+func (z *E4) Add(x, y *E4) *E4 {
+	z.B0.Add(&x.B0, &y.B0)
+	z.B1.Add(&x.B1, &y.B1)
+	return z
+}
+
+// Sub sets z to x-y and returns z
+func (z *E4) Sub(x, y *E4) *E4 {
+	z.B0.Sub(&x.B0, &y.B0)
+	z.B1.Sub(&x.B1, &y.B1)
+	return z
+}
+
+// Double sets z=2*x and returns z
+func (z *E4) Double(x *E4) *E4 {
+	z.B0.Double(&x.B0)
+	z.B1.Double(&x.B1)
+	return z
+}
+
+// Neg negates an E4 element
+func (z *E4) Neg(x *E4) *E4 {
+	z.B0.Neg(&x.B0)
+	z.B1.Neg(&x.B1)
+	return z
+}
+
+// SetRandom used only in tests
+func (z *E4) SetRandom() (*E4, error) {
+	if _, err := z.B0.SetRandom(); err != nil {
+		return nil, err
+	}
+	if _, err := z.B1.SetRandom(); err != nil {
+		return nil, err
+	}
+	return z, nil
+}
+
+// IsZero returns true if z is zero, false otherwise
+func (z *E4) IsZero() bool {
+	return z.B0.IsZero() && z.B1.IsZero()
+}
+
+// IsOne returns true if z is one, false otherwise
+func (z *E4) IsOne() bool {
+	return z.B0.IsOne() && z.B1.IsZero()
+}
+
+// MulByNonResidue mul x by (0,1)
+func (z *E4) MulByNonResidue(x *E4) *E4 {
+	z.B1, z.B0 = x.B0, x.B1
+	z.B0.MulByNonResidue(&z.B0)
+	return z
+}
+
+// Mul sets z=x*y in E4 and returns z
+func (z *E4) Mul(x, y *E4) *E4 {
+	var a, b, c E2
+	a.Add(&x.B0, &x.B1)
+	b.Add(&y.B0, &y.B1)
+	a.Mul(&a, &b)
+	b.Mul(&x.B0, &y.B0)
+	c.Mul(&x.B1, &y.B1)
+	z.B1.Sub(&a, &b).Sub(&z.B1, &c)
+	z.B0.MulByNonResidue(&c).Add(&z.B0, &b)
+	return z
+}
+
+// Square sets z=x*x in E4 and returns z
+func (z *E4) Square(x *E4) *E4 {
+
+	//Algorithm 22 from https://eprint.iacr.org/2010/354.pdf
+	var c0, c2, c3 E2
+	c0.Sub(&x.B0, &x.B1)
+	c3.MulByNonResidue(&x.B1).Sub(&x.B0, &c3)
+	c2.Mul(&x.B0, &x.B1)
+	c0.Mul(&c0, &c3).Add(&c0, &c2)
+	z.B1.Double(&c2)
+	c2.MulByNonResidue(&c2)
+	z.B0.Add(&c0, &c2)
+
+	return z
+}
+
+// Inverse sets z to the inverse of x in E4 and returns z
+//
+// if x == 0, sets and returns z = x
+func (z *E4) Inverse(x *E4) *E4 {
+	// Algorithm 23 from https://eprint.iacr.org/2010/354.pdf
+
+	var t0, t1, tmp E2
+	t0.Square(&x.B0)
+	t1.Square(&x.B1)
+	tmp.MulByNonResidue(&t1)
+	t0.Sub(&t0, &tmp)
+	t1.Inverse(&t0)
+	z.B0.Mul(&x.B0, &t1)
+	z.B1.Mul(&x.B1, &t1).Neg(&z.B1)
+
+	return z
+}
+
+// Exp sets z=xᵏ (mod q⁴) and returns it
+func (z *E4) Exp(x E4, k *big.Int) *E4 {
+	if k.IsUint64() && k.Uint64() == 0 {
+		return z.SetOne()
+	}
+
+	e := k
+	if k.Sign() == -1 {
+		// negative k, we invert
+		// if k < 0: xᵏ (mod q⁴) == (x⁻¹)ᵏ (mod q⁴)
+		x.Inverse(&x)
+
+		// we negate k in a temp big.Int since
+		// Int.Bit(_) of k and -k is different
+		e = bigIntPool.Get().(*big.Int)
+		defer bigIntPool.Put(e)
+		e.Neg(k)
+	}
+
+	z.SetOne()
+	b := e.Bytes()
+	for i := 0; i < len(b); i++ {
+		w := b[i]
+		for j := 0; j < 8; j++ {
+			z.Square(z)
+			if (w & (0b10000000 >> j)) != 0 {
+				z.Mul(z, &x)
+			}
+		}
+	}
+
+	return z
+}
+
+// Conjugate sets z to x conjugated and returns z
+func (z *E4) Conjugate(x *E4) *E4 {
+	z.B0 = x.B0
+	z.B1.Neg(&x.B1)
+	return z
+}
+
+func (z *E4) Halve() {
+
+	z.B0.A0.Halve()
+	z.B0.A1.Halve()
+	z.B1.A0.Halve()
+	z.B1.A1.Halve()
+}
+
+// norm sets x to the norm of z
+func (z *E4) norm(x *E2) {
+	var tmp E2
+	tmp.Square(&z.B1).MulByNonResidue(&tmp)
+	x.Square(&z.B0).Sub(x, &tmp)
+}
+
+// Legendre returns the Legendre symbol of z
+func (z *E4) Legendre() int {
+	var n E2
+	z.norm(&n)
+	return n.Legendre()
+}
+
+// Sqrt sets z to the square root of and returns z
+// The function does not test whether the square root
+// exists or not, it's up to the caller to call
+// Legendre beforehand.
+// cf https://eprint.iacr.org/2012/685.pdf (algo 10)
+func (z *E4) Sqrt(x *E4) *E4 {
+
+	// precomputation
+	var b, c, d, e, f, x0, _g E4
+	var _b, o E2
+
+	// c must be a non square (p = 1 mod 4)
+	c.B1.SetOne()
+
+	q := fr.Modulus()
+	var exp, one big.Int
+	one.SetUint64(1)
+	exp.Mul(q, q).Sub(&exp, &one).Rsh(&exp, 1)
+	d.Exp(c, &exp)
+	e.Mul(&d, &c).Inverse(&e)
+	f.Mul(&d, &c).Square(&f)
+
+	// computation
+	exp.Rsh(&exp, 1)
+	b.Exp(*x, &exp)
+	b.norm(&_b)
+	o.SetOne()
+	if _b.Equal(&o) {
+		x0.Square(&b).Mul(&x0, x)
+		_b.Set(&x0.B0).Sqrt(&_b)
+		_g.B0.Set(&_b)
+		z.Conjugate(&b).Mul(z, &_g)
+		return z
+	}
+	x0.Square(&b).Mul(&x0, x).Mul(&x0, &f)
+	_b.Set(&x0.B0).Sqrt(&_b)
+	_g.B0.Set(&_b)
+	z.Conjugate(&b).Mul(z, &_g).Mul(z, &e)
+
+	return z
+}
+
+// BatchInvertE4 returns a new slice with every element in a inverted.
+// It uses Montgomery batch inversion trick.
+//
+// if a[i] == 0, returns result[i] = a[i]
+func BatchInvertE4(a []E4) []E4 {
+	res := make([]E4, len(a))
+	if len(a) == 0 {
+		return res
+	}
+
+	zeroes := make([]bool, len(a))
+	var accumulator E4
+	accumulator.SetOne()
+
+	for i := 0; i < len(a); i++ {
+		if a[i].IsZero() {
+			zeroes[i] = true
+			continue
+		}
+		res[i].Set(&accumulator)
+		accumulator.Mul(&accumulator, &a[i])
+	}
+
+	accumulator.Inverse(&accumulator)
+
+	for i := len(a) - 1; i >= 0; i-- {
+		if zeroes[i] {
+			continue
+		}
+		res[i].Mul(&res[i], &accumulator)
+		accumulator.Mul(&accumulator, &a[i])
+	}
+
+	return res
+}
+
+// Div divides an element in E4 by an element in E4
+func (z *E4) Div(x *E4, y *E4) *E4 {
+	var r E4
+	r.Inverse(y).Mul(x, &r)
+	return z.Set(&r)
+}
diff --git a/field/koalabear/extensions/e4_test.go b/field/koalabear/extensions/e4_test.go
new file mode 100644
index 0000000000..96f80cc099
--- /dev/null
+++ b/field/koalabear/extensions/e4_test.go
@@ -0,0 +1,300 @@
+package extensions
+
+import (
+	"testing"
+
+	"github.com/leanovate/gopter"
+	"github.com/leanovate/gopter/prop"
+)
+
+// ------------------------------------------------------------
+// tests
+
+func TestE4ReceiverIsOperand(t *testing.T) {
+
+	parameters := gopter.DefaultTestParameters()
+	parameters.MinSuccessfulTests = 100
+
+	properties := gopter.NewProperties(parameters)
+
+	genA := GenE4()
+	genB := GenE4()
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (addition) should output the same result", prop.ForAll(
+		func(a, b *E4) bool {
+			var c, d E4
+			d.Set(a)
+			c.Add(a, b)
+			a.Add(a, b)
+			b.Add(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (sub) should output the same result", prop.ForAll(
+		func(a, b *E4) bool {
+			var c, d E4
+			d.Set(a)
+			c.Sub(a, b)
+			a.Sub(a, b)
+			b.Sub(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (mul) should output the same result", prop.ForAll(
+		func(a, b *E4) bool {
+			var c, d E4
+			d.Set(a)
+			c.Mul(a, b)
+			a.Mul(a, b)
+			b.Mul(&d, b)
+			return a.Equal(b) && a.Equal(&c) && b.Equal(&c)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (square) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Square(a)
+			a.Square(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (double) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Double(a)
+			a.Double(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (mul by non residue) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.MulByNonResidue(a)
+			a.MulByNonResidue(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (Inverse) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Inverse(a)
+			a.Inverse(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (Conjugate) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Conjugate(a)
+			a.Conjugate(a)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Having the receiver as operand (Sqrt) should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b, c, d, s E4
+
+			s.Square(a)
+			a.Set(&s)
+			b.Set(&s)
+
+			a.Sqrt(a)
+			b.Sqrt(&b)
+
+			c.Square(a)
+			d.Square(&b)
+			return c.Equal(&d)
+		},
+		genA,
+	))
+
+	properties.TestingRun(t, gopter.ConsoleReporter(false))
+}
+
+func TestE4Ops(t *testing.T) {
+
+	parameters := gopter.DefaultTestParameters()
+	parameters.MinSuccessfulTests = 100
+
+	properties := gopter.NewProperties(parameters)
+
+	genA := GenE4()
+	genB := GenE4()
+
+	properties.Property("[KOALABEAR] sub & add should leave an element invariant", prop.ForAll(
+		func(a, b *E4) bool {
+			var c E4
+			c.Set(a)
+			c.Add(&c, b).Sub(&c, b)
+			return c.Equal(a)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[KOALABEAR] mul & inverse should leave an element invariant", prop.ForAll(
+		func(a, b *E4) bool {
+			var c, d E4
+			d.Inverse(b)
+			c.Set(a)
+			c.Mul(&c, b).Mul(&c, &d)
+			return c.Equal(a)
+		},
+		genA,
+		genB,
+	))
+
+	properties.Property("[KOALABEAR] BatchInvertE4 should output the same result as Inverse", prop.ForAll(
+		func(a, b, c *E4) bool {
+
+			batch := BatchInvertE4([]E4{*a, *b, *c})
+			a.Inverse(a)
+			b.Inverse(b)
+			c.Inverse(c)
+			return a.Equal(&batch[0]) && b.Equal(&batch[1]) && c.Equal(&batch[2])
+		},
+		genA,
+		genA,
+		genB,
+	))
+
+	properties.Property("[KOALABEAR] inverse twice should leave an element invariant", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Inverse(a).Inverse(&b)
+			return a.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] square and mul should output the same result", prop.ForAll(
+		func(a *E4) bool {
+			var b, c E4
+			b.Mul(a, a)
+			c.Square(a)
+			return b.Equal(&c)
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] Legendre on square should output 1", prop.ForAll(
+		func(a *E4) bool {
+			var b E4
+			b.Square(a)
+			c := b.Legendre()
+			return c == 1
+		},
+		genA,
+	))
+
+	properties.Property("[KOALABEAR] square(sqrt) should leave an element invariant", prop.ForAll(
+		func(a *E4) bool {
+			var b, c, d, e E4
+			b.Square(a)
+			c.Sqrt(&b)
+			d.Square(&c)
+			e.Neg(a)
+			return (c.Equal(a) || c.Equal(&e)) && d.Equal(&b)
+		},
+		genA,
+	))
+
+	properties.TestingRun(t, gopter.ConsoleReporter(false))
+}
+
+// ------------------------------------------------------------
+// benches
+
+func BenchmarkE4Add(b *testing.B) {
+	var a, c E4
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Add(&a, &c)
+	}
+}
+
+func BenchmarkE4Sub(b *testing.B) {
+	var a, c E4
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Sub(&a, &c)
+	}
+}
+
+func BenchmarkE4Mul(b *testing.B) {
+	var a, c E4
+	_, _ = a.SetRandom()
+	_, _ = c.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Mul(&a, &c)
+	}
+}
+
+func BenchmarkE4Square(b *testing.B) {
+	var a E4
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Square(&a)
+	}
+}
+
+func BenchmarkE4Sqrt(b *testing.B) {
+	var a E4
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Sqrt(&a)
+	}
+}
+
+func BenchmarkE4Inverse(b *testing.B) {
+	var a E4
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Inverse(&a)
+	}
+}
+
+func BenchmarkE4MulNonRes(b *testing.B) {
+	var a E4
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.MulByNonResidue(&a)
+	}
+}
+
+func BenchmarkE4Conjugate(b *testing.B) {
+	var a E4
+	_, _ = a.SetRandom()
+	b.ResetTimer()
+	for i := 0; i < b.N; i++ {
+		a.Conjugate(&a)
+	}
+}
diff --git a/field/koalabear/extensions/utils.go b/field/koalabear/extensions/utils.go
new file mode 100644
index 0000000000..998332f1e1
--- /dev/null
+++ b/field/koalabear/extensions/utils.go
@@ -0,0 +1,48 @@
+package extensions
+
+import (
+	"math/big"
+	"sync"
+
+	fr "github.com/consensys/gnark-crypto/field/koalabear"
+	"github.com/leanovate/gopter"
+)
+
+var bigIntPool = sync.Pool{
+	New: func() interface{} {
+		return new(big.Int)
+	},
+}
+
+// Fr generates an Fr element
+func GenFr() gopter.Gen {
+	return func(genParams *gopter.GenParameters) *gopter.GenResult {
+		var elmt fr.Element
+
+		if _, err := elmt.SetRandom(); err != nil {
+			panic(err)
+		}
+		genResult := gopter.NewGenResult(elmt, gopter.NoShrinker)
+		return genResult
+	}
+}
+
+// E2 generates an E2 elmt
+func GenE2() gopter.Gen {
+	return gopter.CombineGens(
+		GenFr(),
+		GenFr(),
+	).Map(func(values []interface{}) *E2 {
+		return &E2{A0: values[0].(fr.Element), A1: values[1].(fr.Element)}
+	})
+}
+
+// E4 generates an E4 elmt
+func GenE4() gopter.Gen {
+	return gopter.CombineGens(
+		GenE2(),
+		GenE2(),
+	).Map(func(values []interface{}) *E4 {
+		return &E4{B0: *values[0].(*E2), B1: *values[1].(*E2)}
+	})
+}