Hi folks, Can I just get a sanity check on this bit of code to check I'm not doing something insane (I probably am) or stupid (equally likely). The problem is a basic diffusion equation where a pile of sediment is eroded and transport by diffusion only. The un-erodable topographic surface (land) + sediment (phi) make up the land surface and the slope over which the sediment is moved. The attached figure might show it more clearly. phi is therefore moved according to the 2nd order (phi+land). However, I've also had to "turn off" diffusion by setting the diffusion coefficient (nu) to zero where there is no sediment or volume is definitely not conserved. The equation is set up using: tiny = 1e-10 surface.interpolate((((land+phi) + land)+abs((land+phi)-land)) / 2) thickness.interpolate(surface-land) limit.interpolate((surface - land)/(surface-land+tiny)) F = (inner((phi - phi_)/timestep, v) + limit*nu*inner(grad(phi+land), grad(v)))*dx And then simply solved: solve(F == 0, phi) in the time loop, updating limit and surface as we go. Hopefully, the above is understandable, but if anyone has any opinions on my (in)sanity, then I'm happy to hear them! Cheers, J. -- I'm currently in Australia and replies will not be in standard UK business hours. Please do not feel you have to respond outside your normal work hours. Dr Jon Hill Senior Lecturer in Physical Geography Department of Environment and Geography University of York M: +44(0)7748254812 P: +44(0)1904 324480 Twitter: @jonxhill Web: https://jonxhill.wordpress.com/ Web: https://envmodellinggroup.github.io/
