CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > General Forums > Main CFD Forum

Solving u_t = Au_x +Bu_y

Register Blogs Community New Posts Updated Threads Search

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   September 8, 2011, 09:51
Default Solving u_t = A(u_x +u_y) in MATLAB
  #1
New Member
 
Vincent
Join Date: Jul 2011
Posts: 29
Rep Power: 14
VincentD is on a distinguished road
I would like to solve the following equation:
u_t =A(u_x + u_y)
Where A is a matrix with size equal to u. I use MATLAB to do the computation.
I want the scheme to be explicit in time and implicit in space.

Before I adress this equation I will first make it a step more simple.
Lets say we want to solve
u_t = a(u_x +u_y)
Where a is a scalar.

Rewriting this equation gives:

u*-dt a(u*_x +u*_y) = u^n

We can write this equation as a matrix product in the following way:

(I+dt a(D2D))u* =u^n

Where I is the identity matrix and D2D is the matrix approximating -d/dx -d/dy. We would this D2D matrix in the following way.

We define the difference matrix D (approximating -d/dn, 1D)

spdiags([1 1 -1 ;ones(n-2,1)*[1 0 -1]; 1 1 -1],-1:1,n,n)'/(2*h); (Neumann BCs)

And we can create D2D by:

kron(speye(ny),D2D(nx-1))+kron(D2D(ny),speye(nx-1))

The full matrix I+dt a(D2D) is than

kron(speye(nx-1),speye(ny))+dt*a*kron(speye(ny),D(nx-1,hx,2,1))+kron(D(ny,hy,3,3),speye(nx-1))) (don't mind the arguments in D, first one denotes size)

Alright, well this works. But now we make it a little bit more difficult, namely we change a to be a matrix instead of a scalar

u_t =A(u_x + u_y)

The matrix equation becomes
(I+dt A(D2D))u* =u^n

Now we take that nx=ny=10
Now A is a matrix of 9x10 and D2D has size 90x90 and same for I (90x90).

I need to blow matrix A up in such a way that it's also square with a size (90x90), however if I agian use the kron command

kron(speye(10),A) than the size is 100x90. If matrix A would be square than a multiplication like this would certainly work, but alas it is not.

If anyone has a clue on how to do this I would love to hear it.

Regards,

Vincent

Last edited by VincentD; September 9, 2011 at 05:02.
VincentD is offline   Reply With Quote

Reply


Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Moving mesh Niklas Wikstrom (Wikstrom) OpenFOAM Running, Solving & CFD 122 June 15, 2014 06:20
Upgraded from Karmic Koala 9.10 to Lucid Lynx10.04.3 bookie56 OpenFOAM Installation 8 August 13, 2011 04:03
Orifice Plate with a fully developed flow - Problems with convergence jonmec OpenFOAM Running, Solving & CFD 3 July 28, 2011 05:24
Differences between serial and parallel runs carsten OpenFOAM Bugs 11 September 12, 2008 11:16
Could anybody help me see this error and give help liugx212 OpenFOAM Running, Solving & CFD 3 January 4, 2006 18:07


All times are GMT -4. The time now is 01:07.