diff --git a/geom.go b/geom.go
index 387592c70..f199c59ce 100644
--- a/geom.go
+++ b/geom.go
@@ -26,12 +26,12 @@ const GeoMDim = 3
 //
 // The initial value is identity.
 type GeoM struct {
-	a_1 float64 // The actual 'a' value minus 1
-	b   float64
-	c   float64
-	d_1 float64 // The actual 'd' value minus 1
-	tx  float64
-	ty  float64
+	a_1 float32 // The actual 'a' value minus 1
+	b   float32
+	c   float32
+	d_1 float32 // The actual 'd' value minus 1
+	tx  float32
+	ty  float32
 }
 
 // String returns a string representation of GeoM.
@@ -52,7 +52,12 @@ func (g *GeoM) Reset() {
 // Apply pre-multiplies a vector (x, y, 1) by the matrix.
 // In other words, Apply calculates GeoM * (x, y, 1)^T.
 // The return value is x and y values of the result vector.
-func (g *GeoM) Apply(x, y float64) (x2, y2 float64) {
+func (g *GeoM) Apply(x, y float64) (float64, float64) {
+	x2, y2 := g.apply32(float32(x), float32(y))
+	return float64(x2), float64(y2)
+}
+
+func (g *GeoM) apply32(x, y float32) (x2, y2 float32) {
 	return (g.a_1+1)*x + g.b*y + g.tx, g.c*x + (g.d_1+1)*y + g.ty
 }
 
@@ -60,17 +65,17 @@ func (g *GeoM) Apply(x, y float64) (x2, y2 float64) {
 func (g *GeoM) Element(i, j int) float64 {
 	switch {
 	case i == 0 && j == 0:
-		return g.a_1 + 1
+		return float64(g.a_1) + 1
 	case i == 0 && j == 1:
-		return g.b
+		return float64(g.b)
 	case i == 0 && j == 2:
-		return g.tx
+		return float64(g.tx)
 	case i == 1 && j == 0:
-		return g.c
+		return float64(g.c)
 	case i == 1 && j == 1:
-		return g.d_1 + 1
+		return float64(g.d_1) + 1
 	case i == 1 && j == 2:
-		return g.ty
+		return float64(g.ty)
 	default:
 		panic("ebiten: i or j is out of index")
 	}
@@ -107,31 +112,32 @@ func (g *GeoM) Add(other GeoM) {
 
 // Scale scales the matrix by (x, y).
 func (g *GeoM) Scale(x, y float64) {
-	a := (g.a_1 + 1) * x
-	b := g.b * x
-	tx := g.tx * x
-	c := g.c * y
-	d := (g.d_1 + 1) * y
-	ty := g.ty * y
+	a := (float64(g.a_1) + 1) * x
+	b := float64(g.b) * x
+	tx := float64(g.tx) * x
+	c := float64(g.c) * y
+	d := (float64(g.d_1) + 1) * y
+	ty := float64(g.ty) * y
 
-	g.a_1 = a - 1
-	g.b = b
-	g.c = c
-	g.d_1 = d - 1
-	g.tx = tx
-	g.ty = ty
+	g.a_1 = float32(a) - 1
+	g.b = float32(b)
+	g.c = float32(c)
+	g.d_1 = float32(d) - 1
+	g.tx = float32(tx)
+	g.ty = float32(ty)
 }
 
 // Translate translates the matrix by (tx, ty).
 func (g *GeoM) Translate(tx, ty float64) {
-	g.tx += tx
-	g.ty += ty
+	g.tx += float32(tx)
+	g.ty += float32(ty)
 }
 
 // Rotate rotates the matrix by theta.
 // The unit is radian.
 func (g *GeoM) Rotate(theta float64) {
-	sin, cos := math.Sincos(theta)
+	sin64, cos64 := math.Sincos(theta)
+	sin, cos := float32(sin64), float32(cos64)
 
 	a := cos*(g.a_1+1) - sin*g.c
 	b := cos*g.b - sin*(g.d_1+1)
@@ -148,7 +154,7 @@ func (g *GeoM) Rotate(theta float64) {
 	g.ty = ty
 }
 
-func (g *GeoM) det() float64 {
+func (g *GeoM) det() float32 {
 	return (g.a_1+1)*(g.d_1+1) - g.b*g.c
 }
 
@@ -183,19 +189,20 @@ func (g *GeoM) Invert() {
 
 // SetElement sets an element at (i, j).
 func (g *GeoM) SetElement(i, j int, element float64) {
+	e := float32(element)
 	switch {
 	case i == 0 && j == 0:
-		g.a_1 = element - 1
+		g.a_1 = e - 1
 	case i == 0 && j == 1:
-		g.b = element
+		g.b = e
 	case i == 0 && j == 2:
-		g.tx = element
+		g.tx = e
 	case i == 1 && j == 0:
-		g.c = element
+		g.c = e
 	case i == 1 && j == 1:
-		g.d_1 = element - 1
+		g.d_1 = e - 1
 	case i == 1 && j == 2:
-		g.ty = element
+		g.ty = e
 	default:
 		panic("ebiten: i or j is out of index")
 	}
diff --git a/geom_test.go b/geom_test.go
index bdb289d8a..32d666378 100644
--- a/geom_test.go
+++ b/geom_test.go
@@ -258,7 +258,7 @@ func TestGeoMIsInvert(t *testing.T) {
 		},
 	}
 
-	const delta = 0.00001
+	const delta = 0.001
 
 	for _, c := range cases {
 		if c.GeoM.IsInvertible() != c.Invertible {
diff --git a/image.go b/image.go
index 2a8111733..5b3df3e2d 100644
--- a/image.go
+++ b/image.go
@@ -126,6 +126,14 @@ func (i *Image) fill(r, g, b, a uint8) {
 	_ = i.DrawImage(emptyImage, op)
 }
 
+type geoM32 struct {
+	inner *GeoM
+}
+
+func (g geoM32) Apply(x, y float32) (x2, y2 float32) {
+	return g.inner.apply32(x, y)
+}
+
 // DrawImage draws the given image on the image i.
 //
 // DrawImage accepts the options. For details, see the document of DrawImageOptions.
@@ -211,7 +219,7 @@ func (i *Image) DrawImage(img *Image, options *DrawImageOptions) error {
 			sy1 = r.Max.Y
 		}
 	}
-	geom := &options.GeoM
+	geom := geoM32{&options.GeoM}
 	if sx0 < 0 || sy0 < 0 {
 		dx := 0.0
 		dy := 0.0
@@ -223,9 +231,10 @@ func (i *Image) DrawImage(img *Image, options *DrawImageOptions) error {
 			dy = -float64(sy0)
 			sy0 = 0
 		}
-		geom = &GeoM{}
-		geom.Translate(dx, dy)
-		geom.Concat(options.GeoM)
+		g := &GeoM{}
+		g.Translate(dx, dy)
+		g.Concat(options.GeoM)
+		geom = geoM32{g}
 	}
 
 	mode := opengl.CompositeMode(options.CompositeMode)
diff --git a/internal/graphicsutil/vertices.go b/internal/graphicsutil/vertices.go
index 7a2e9f723..d5414d987 100644
--- a/internal/graphicsutil/vertices.go
+++ b/internal/graphicsutil/vertices.go
@@ -44,7 +44,7 @@ func (v *verticesBackend) sliceForOneQuad() []float32 {
 }
 
 type GeoM interface {
-	Apply(x, y float64) (x2, y2 float64)
+	Apply(x, y float32) (x2, y2 float32)
 }
 
 func QuadVertices(width, height int, sx0, sy0, sx1, sy1 int, geom GeoM) []float32 {
@@ -57,8 +57,8 @@ func QuadVertices(width, height int, sx0, sy0, sx1, sy1 int, geom GeoM) []float3
 
 	vs := theVerticesBackend.sliceForOneQuad()
 
-	x0, y0 := 0.0, 0.0
-	x1, y1 := float64(sx1-sx0), float64(sy1-sy0)
+	x0, y0 := float32(0.0), float32(0.0)
+	x1, y1 := float32(sx1-sx0), float32(sy1-sy0)
 
 	// it really feels like we should be able to cache this computation
 	// but it may not matter.
@@ -73,11 +73,15 @@ func QuadVertices(width, height int, sx0, sy0, sx1, sy1 int, geom GeoM) []float3
 	wf := float32(w)
 	hf := float32(h)
 	u0, v0, u1, v1 := float32(sx0)/wf, float32(sy0)/hf, float32(sx1)/wf, float32(sy1)/hf
+	quadVerticesImpl(vs, wf, hf, u0, v0, u1, v1, x0, y0, x1, y1, geom)
+	return vs
+}
 
+func quadVerticesImpl(vs []float32, wf, hf, u0, v0, u1, v1, x0, y0, x1, y1 float32, geom GeoM) {
 	x, y := geom.Apply(x0, y0)
 	// Vertex coordinates
-	vs[0] = float32(x)
-	vs[1] = float32(y)
+	vs[0] = x
+	vs[1] = y
 
 	// Texture coordinates: first 2 values indicates the actual coodinate, and
 	// the second indicates diagonally opposite coodinates.
@@ -89,28 +93,26 @@ func QuadVertices(width, height int, sx0, sy0, sx1, sy1 int, geom GeoM) []float3
 
 	// and the same for the other three coordinates
 	x, y = geom.Apply(x1, y0)
-	vs[6] = float32(x)
-	vs[7] = float32(y)
+	vs[6] = x
+	vs[7] = y
 	vs[8] = u1
 	vs[9] = v0
 	vs[10] = u0
 	vs[11] = v1
 
 	x, y = geom.Apply(x0, y1)
-	vs[12] = float32(x)
-	vs[13] = float32(y)
+	vs[12] = x
+	vs[13] = y
 	vs[14] = u0
 	vs[15] = v1
 	vs[16] = u1
 	vs[17] = v0
 
 	x, y = geom.Apply(x1, y1)
-	vs[18] = float32(x)
-	vs[19] = float32(y)
+	vs[18] = x
+	vs[19] = y
 	vs[20] = u1
 	vs[21] = v1
 	vs[22] = u0
 	vs[23] = v0
-
-	return vs
 }
diff --git a/internal/restorable/image.go b/internal/restorable/image.go
index 012928e53..762b40cdc 100644
--- a/internal/restorable/image.go
+++ b/internal/restorable/image.go
@@ -123,13 +123,13 @@ var (
 )
 
 type geoM struct {
-	scaleX float64
-	scaleY float64
-	tx     float64
-	ty     float64
+	scaleX float32
+	scaleY float32
+	tx     float32
+	ty     float32
 }
 
-func (g *geoM) Apply(x, y float64) (float64, float64) {
+func (g *geoM) Apply(x, y float32) (float32, float32) {
 	return x*g.scaleX + g.tx, y*g.scaleY + g.ty
 }
 
@@ -159,10 +159,10 @@ func (i *Image) ReplacePixels(pixels []byte, x, y, width, height int) {
 		w, h := dummyImage.Size()
 		colorm := (*affine.ColorM)(nil).Scale(0, 0, 0, 0)
 		vs := graphicsutil.QuadVertices(w, h, 0, 0, w, h, &geoM{
-			scaleX: float64(width) / float64(w),
-			scaleY: float64(height) / float64(h),
-			tx:     float64(x),
-			ty:     float64(y),
+			scaleX: float32(width) / float32(w),
+			scaleY: float32(height) / float32(h),
+			tx:     float32(x),
+			ty:     float32(y),
 		})
 		i.image.DrawImage(dummyImage.image, vs, quadIndices, colorm, opengl.CompositeModeCopy, graphics.FilterNearest)
 	}
diff --git a/internal/restorable/images_test.go b/internal/restorable/images_test.go
index a82140482..cdd037a32 100644
--- a/internal/restorable/images_test.go
+++ b/internal/restorable/images_test.go
@@ -107,7 +107,7 @@ var (
 
 type idGeoM struct{}
 
-func (idGeoM) Apply(x, y float64) (x2, y2 float64) {
+func (idGeoM) Apply(x, y float32) (x2, y2 float32) {
 	return x, y
 }
 
@@ -263,11 +263,11 @@ func TestRestoreOverrideSource(t *testing.T) {
 }
 
 type geoM struct {
-	tx float64
-	ty float64
+	tx float32
+	ty float32
 }
 
-func (g *geoM) Apply(x, y float64) (x2, y2 float64) {
+func (g *geoM) Apply(x, y float32) (x2, y2 float32) {
 	return x + g.tx, y + g.ty
 }
 
diff --git a/internal/shareable/shareable.go b/internal/shareable/shareable.go
index 143a83bf6..fd83b1947 100644
--- a/internal/shareable/shareable.go
+++ b/internal/shareable/shareable.go
@@ -100,7 +100,7 @@ type Image struct {
 
 type idGeoM struct{}
 
-func (idGeoM) Apply(x, y float64) (x2, y2 float64) {
+func (idGeoM) Apply(x, y float32) (x2, y2 float32) {
 	return x, y
 }
 
diff --git a/internal/shareable/shareable_test.go b/internal/shareable/shareable_test.go
index b22d7fa1d..b2a929959 100644
--- a/internal/shareable/shareable_test.go
+++ b/internal/shareable/shareable_test.go
@@ -46,11 +46,11 @@ func TestMain(m *testing.M) {
 const bigSize = 2049
 
 type geoM struct {
-	tx float64
-	ty float64
+	tx float32
+	ty float32
 }
 
-func (g *geoM) Apply(x, y float64) (x2, y2 float64) {
+func (g *geoM) Apply(x, y float32) (x2, y2 float32) {
 	return x + g.tx, y + g.ty
 }