Looking at topdown 2D car physics code and tutorials, a lot of them cite Marco Monster’s article (archived here, a port of the code here). Here’s a quick Codea port I made, just the car/wheel physics (not the engine, gears etc), using modelMatrix
to handle the rotations. The basic principle is that you translate into the local space of the car, then you can separate the lateral and longitudinal forces acting on the body and the wheels.
I don’t understand how all the variables in the simulation interact, but you can make the car drift quite easily.
You could also implement this in Box2D, by cancelling out most of the lateral forces on the wheels, eg top example here.
Tilt to steer, touch right side to accelerate, left side to brake/ reverse.
Issues:

Friction doesn’t seem to bring car to a complete halt, there seems to be lingering velocity.

Steering doesn’t work when vehicle is going backwards
# Main
 2D Car Physics
supportedOrientations(LANDSCAPE_ANY)
displayMode(OVERLAY)
function setup()
centre=vec2(WIDTH,HEIGHT)/2
parameter.watch("you.pos")
parameter.watch("you.angleVel")
parameter.watch("you.vel")
parameter.watch("you.steerAngle")
you = Vehicle{pos=vec2(WIDTH,HEIGHT)/2, size=vec2(50,30), mass=700, inertia=700, steering = Vehicle.tiltSteering, gas = Vehicle.touchGas}
end
function draw()
background(40, 40, 50)
you:update()
end
function touched(t)
if t.state == ENDED then
touching = false
else
if t.x<centre.x then
touching = 1
else
touching = 2
end
end
end
# Body
Body = class()
function Body:init(t)
self.size=t.size
local w, h = t.size.x/2, t.size.y/2
self.radius = vec2(w,h)
self.points={vec2(w,h), vec2(w,h), vec2(w,h), vec2(w,h)}
 linear properties
self.pos = t.pos
self.vel = t.vel or vec2(0,0)
self.forces = vec2(0,0) accumulated forces. Zeroed every frame
self.mass = t.mass or 1 kg
angular properties
self.angle = t.angle or 0
self.angleVel = 0 in rad/s
self.torque = 0 accumulated torque. zeroed every frame
self.inertia = t.inertia or 150 kg.m
total transform
self.matrix=matrix()
end
function Body:constrain()
self.vel = outOfBounds(self.vel, self.posself.radius, self.pos+self.radius)
end
function Body:integrate()
local dt = DeltaTime
integrate physics
linear
local accel = self.forces /self.mass
self.vel = self.vel + dt * accel
angular
local angAcc = self.torque /self.inertia
self.angleVel = self.angleVel + dt * angAcc
constrain
 self:constrain()
update
self.pos = self.pos + dt *self.vel
self.angle = self.angle + dt * self.angleVel
reset accumulators
self.forces = vec2(0,0)
self.torque = 0
end
function Body:getWorldPoint(v)
return self.matrix * v
end
function Body:getLocalPoint(v)
local inv = self.matrix:inverse()
return inv * v
end
function Body:applyForce(force, worldPoint)
self.forces = self.forces + force
[[
if worldPoint then
self:applyTorque( math.rad( worldPointself.pos:cross(force)) ) world offset (torque arm) crossed with world force. Why does this have to be converted to radians though?
end
]]
end
function Body:applyTorque(ang)
self.torque = self.torque + ang
end
# Vehicle
Vehicle = class(Body) based on Marco Monster's model
X Forwards.
added weight distribution formula
issues: steering doesnt work in reverse
friction doesnt quite bring car to complete stop
local drag = 20 4.257 air resistance, proportional to vel squared
local resistance = drag * 30 rolling resistance, proportional to vel
local cornerStiffRear = 5.2
local cornerStiffFront = 5
local maxGrip= 4 increase this for a tighter turning circle, less sliding?
local engineAccel = 180
local engineBreak = engineAccel * 10
local sensitivity = 1 steering sensitivity
function Vehicle:init(t)
Body.init(self, t)
self.wheelRadius = 10
self.length = t.size.x  (self.wheelRadius*2) wheelbase
self.halfLength = self.length/2
weight distribution variables
local height=2 vertical distance from centre of wheels to centre of geometry (keep low, otherwise vehicle really slides on breaking)
self.heightRatio = 2/ self.length  height as a proportion if wheelbase, used to calculate shifting weight distribution
self.weightDiff = 0
self.accelRate = 0  this is stored to calculate distribution of weight between axles
inputs
self.steering = t.steering
self.steerTarget=0
self.gas = t.gas
self.revs = 0
self.maxRevs = 64000 completely madeup number
end
function Vehicle:update()
 inputs
self:gas()
self:steering()
physics
self:doPhysics()
self:constrain()
self:integrate()
self:draw()
end
function Vehicle:doPhysics()
transform world velocity to local 2D space
local vel = vecRotMat(self.vel, self.matrix:inverse()) 
 local vel = vec2(cos * self.vel.y + sin * self.vel.x, sin * self.vel.y + cos * self.vel.x)
 local vel=returnVec2(self:getLocalPoint(self.vel))
lateral force on wheels. r * vel in radians/second
local yawSpeed = self.halfLength * self.angleVel
rotAngle = math.atan(yawSpeed,vel.x)
sideslip angle of car (beta). ie angle between orientation and velocity
sideSlip = math.atan(vel.y, vel.x)
calculate slip angles for front and rear wheels (alpha)
local slipAngleFront = sideSlip + rotAngle  self.steerAngle
local slipAngleRear = sideSlip  rotAngle
weight per axle = half car mass times 1gee (9.8 m/s^2)
local weight = self.mass * 4.9
weightDiff = self.heightRatio * self.mass * self.accelRate weight distribution between axles (stored to animate body)
local weightFront = weight  weightDiff
local weightRear = weight + weightDiff
vary this according to acceleration
lateral force on front wheels (corner stiffness * slip angle, clamped to friction circle * load), front slips
local latForceFront = vec2(0,clamp(cornerStiffFront * slipAngleFront, maxGrip, maxGrip)*weightFront)
lateral force on rear wheels
local latForceRear = vec2(0, clamp(cornerStiffRear * slipAngleRear, maxGrip, maxGrip)* weightRear)
longitudinal force on rear wheels = traction. 100 * (self.revs  self.brake * sign(vel.x))
local traction = vec2(100 * self.revs, 0 )
drag and rolling resistance
local friction = vel * (resistance + drag * vel:len())
sum forces
local cornering = (latForceFront):rotate(self.steerAngle) + latForceRear
local force = traction + cornering + friction
APPLY FORCE AND TORQUE
sum torque from lat forces
self:applyTorque(latForceFront.y  latForceRear.y)
convert force back to 3D world space
acceleration = vecRotMat(force, self.matrix)  self:getWorldPoint(accel)
 acceleration = vecAddZ(vec2(cos * accel.y + sin * accel.x, sin * accel.y + cos * accel.x), 0)
self:applyForce(acceleration)
self.accelRate = DeltaTime * force.x / self.mass store acceleration rate to calculate weight distribution next frame
end
function Vehicle:draw()
pushMatrix()
translate(self.pos.x, self.pos.y)
rotate(math.deg(self.angle))
self.matrix = modelMatrix()
sprite("Platformer Art:Block Grass", 0,0,self.size.x, self.size.y)
self:drawWheel(1, self.wheelRadius)
self:drawWheel(2, self.wheelRadius)
self:drawWheel(3, self.wheelRadius, self.steerAngle)
self:drawWheel(4, self.wheelRadius, self.steerAngle)
popMatrix()
end
function Vehicle:drawWheel(n, offset, turn)
local turn = turn or 0
pushMatrix()
translate(self.points[n].x+offset, self.points[n].y)
rotate(math.deg(turn))
sprite("Platformer Art:Block Special", 0,0,20,5)
popMatrix()
end
function Vehicle:tiltSteering()
self.steerAngle = clamp(Gravity.x * sensitivity, 0.40, 0.40)
end
function Vehicle:touchGas()
if touching==2 then
self.revs = math.min(self.revs + engineAccel, self.maxRevs)
elseif touching==1 then
self.revs= math.max(self.revs  engineAccel, self.maxRevs * 0.1)
else
self.revs= math.max(self.revs  engineBreak, 0)
end
end
# Helper
function outOfBounds(vel, a, b, bb, aa) if box ab on heading vel is outside of boundary aabb, returns a corrective vel, plus how much out of bounds the box is
local b=b or a
local aa=aa or vec2(0,0)
local bb=bb or vec2(WIDTH,HEIGHT)
local x,y=0,0
if a.x<aa.x and vel.x<0 then x=a.xaa.x vel.x=vel.x
elseif b.x>bb.x and vel.x>0 then x=b.xbb.x vel.x=vel.x
end
if a.y<aa.y and vel.y<0 then y=a.yaa.y vel.y=vel.y
elseif b.y>bb.y and vel.y>0 then y=b.ybb.y vel.y=vel.y
end
return vel, x, y
end
function AABB(v,aa,bb)
if v.x>aa.x and v.x<bb.x and v.y>aa.y and v.y<bb.y then return true end
end
function sign(x) returns 1 if x<0, else 1
return (x<0 and 1) or 1
end
function clamp(v,low,high)
return math.min(math.max(v, low), high)
end
function vecRotMat(v, m) apply rotation part (top 2x2) of a mat4 to a vec2
return vec2(
m[1]*v.x + m[5]*v.y,
m[2]*v.x + m[6]*v.y)
end