[flotus1]

You have several options for the implementation of a first order Maxwell boundary condition.

The easiest (for finite volumes) would be the way Fluent does it. They obtain the wall-normal gradients for velocity and temperature from finite differences between the values at the wall and a the first cell center.
Then they apply the value calculated from the Maxwell boundary condition to the boundary.
So you need a temporal storage to store the boundary value from the last iteration to use it in the finite differences for the next iteration.

Keep in mind that in doing so, the spatial accuracy of the whole computation is reduced to first order. To prevent this, you could for example evaluate the wall-normal gradient at the wall with a higher order finite difference scheme that uses additional points further away from the wall.

Some kind of under-relaxation might be necessary to keep this kind of BC stable, especially with higher order finite difference schemes for the gradients.