Estou tentando usarNDSolve resolver equações de onda para verificar se é mais fácil e / ou mais rápido usá-lo em vez das minhas características antigas eq. implementação de método.

Estou tendo muita instabilidade que não entendo com o método das características, e como essas são equações simples, eu me pergunto o que está errado ... (espero, não o aspecto físico do problema ...)

ans = Flatten@NDSolve[{
u[t, x]*D[d[t, x], x] + d[t, x]*D[u[t, x], x] + D[d[t, x], t] == 0,
D[d[t, x], x] + u[t, x]/9.8*D[u[t, x], x] + 
 1/9.8*D[u[t, x], t] + 0.0001 u[t, x]*Abs[u[t, x]] == 0,
u[0, x] == 0,
d[0, x] == 3 + x/1000*1,
u[t, 0] == 0,
u[t, 1000] == 0
},
d, {t, 0, 1000}, {x, 0, 1000}, DependentVariables -> {u, d}
]

Animate[Plot[(d /. ans)[t, x], {x, 0, 1000}, 
        PlotRange -> {{0, 1000}, {0, 6}}], {t, 0, 1000}
]

Alguém pode me ajudar

EDITAR

Coloquei oNDSolve solution (após a edição do JxB) com minhas soluções de características, juntas na mesma animação. Eles combinam perto o suficiente, com exceção das oscilações rápidas iniciais. Com o tempo, eles tendem a começar a dessincronizar, mas acredito que isso se deve provavelmente a uma pequena simplificação que temos que admitir ao deduzir as característica

Red:NDsolve; Azul: método das características "manuais";

press F5 (atualize o navegador), para reiniciar a animação det=0.

(escala xx é o número de pontos que usei no meu método "manual", em que cada ponto representa 20 unidades doNDSolve / escala física)

Brincando comNDSolve grid sampling, produz efeitos de oscilação completamente diferentes. Alguém tem ou conhece alguma técnica para garantir uma integração adequada?