# problem with a coupled system

I am very sorry to have to ask this but can you please have a look at my code? The BC's are completely off for this coupled system, nothing like what I would expect it to be.

from cbc.pdesys import *

mesh = UnitSquare(20, 20)

# Change desired items in the problem_parameters dict from cbc.pdesys

problem_

problem = Problem(mesh, problem_parameters)

#solver_

solver_

solver_

coupled = PDESystem([['Pw', 'Sw']], problem, solver_parameters) # Creates FunctionSpace, Functions etc.

class SS2Ph(PDESubSys

def form(self, Pw, v_Pw, Pw_, Sw, Sw_, Sw_1, v_Sw, dt, **kwargs):

FPw = inner(grad(

FSw = inner(grad(

return FPw + FSw

# setup the BCs

def left(x):

return x[0] < DOLFIN_EPS

def right(x):

return x[0] > 1.0-DOLFIN_EPS

def top(x):

return x[1] > 1.0-DOLFIN_EPS

def bottom(x):

return x[1] < DOLFIN_EPS

bcs = [DirichletBC(

# add the pde subsystem

coupled.

class InitialConditio

def __init__(self):

print "initialized initial conditions"

def eval(self, values, x):

values[0] = 10.0

values[1] = 0.5

def value_shape(self):

return (2,)

def update(self):

plot(

plot(

Problem.update = update

problem.q0 = {'PwSw': InitialConditio

problem.prm['T'] = 3.0

problem.

problem.solve()

## Question information

- Language:
- English Edit question

- Status:
- Solved

- For:
- CBC.PDESys Edit question

- Assignee:
- No assignee Edit question

- Solved by:
- Yaakoub El Khamra

- Solved:

- Last query:

- Last reply:

## This question was reopened

Hi Yakoub,

Not quite sure what you mean, it seems to work for me? What's the problem? Think you should use DirichletBC(

Mikael

Den Oct 22, 2012 kl. 6:51 PM skrev Yaakoub El Khamra:

> New question #211982 on CBC.PDESys:

> https:/

>

> I am very sorry to have to ask this but can you please have a look at my code? The BC's are completely off for this coupled system, nothing like what I would expect it to be.

>

> from cbc.pdesys import *

>

> mesh = UnitSquare(20, 20)

>

> # Change desired items in the problem_parameters dict from cbc.pdesys

> problem_

> problem = Problem(mesh, problem_parameters)

>

> #solver_

> solver_

> solver_

> coupled = PDESystem([['Pw', 'Sw']], problem, solver_parameters) # Creates FunctionSpace, Functions etc.

>

>

> class SS2Ph(PDESubSys

> def form(self, Pw, v_Pw, Pw_, Sw, Sw_, Sw_1, v_Sw, dt, **kwargs):

>

> FPw = inner(grad(

> FSw = inner(grad(

> return FPw + FSw

>

>

> # setup the BCs

> def left(x):

> return x[0] < DOLFIN_EPS

> def right(x):

> return x[0] > 1.0-DOLFIN_EPS

> def top(x):

> return x[1] > 1.0-DOLFIN_EPS

> def bottom(x):

> return x[1] < DOLFIN_EPS

>

> bcs = [DirichletBC(

> DirichletBC(

>

>

> # add the pde subsystem

> coupled.

>

> class InitialConditio

> def __init__(self):

> print "initialized initial conditions"

> def eval(self, values, x):

> values[0] = 10.0

> values[1] = 0.5

> def value_shape(self):

> return (2,)

>

> def update(self):

> plot(self.

> plot(self.

>

> Problem.update = update

>

> problem.q0 = {'PwSw': InitialConditio

>

> problem.prm['T'] = 3.0

> problem.

> problem.solve()

>

>

> --

> You received this question notification because you are an answer

> contact for CBC.PDESys.

Yaakoub El Khamra (yelkhamra) said : | #2 |

Well for one thing I would expect the boundary on the left to be 15 (x[0]=0), the right to be 5(x[0]=1). It looks as if the bottom right corner has a value of 5, i.e. as if the BC is rotated by -45 degrees. x[0]=1 and x[1]> 0.5 seems to have a value of 15 (should be 5)

This would be easier if I can attach a screen shot, but I cannot find a way to do it.

Yaakoub El Khamra (yelkhamra) said : | #3 |

Sorry, accidentally hit problem solved.

Den Oct 22, 2012 kl. 7:21 PM skrev Yaakoub El Khamra:

> Question #211982 on CBC.PDESys changed:

> https:/

>

> Status: Answered => Solved

>

> Yaakoub El Khamra confirmed that the question is solved:

>

> Well for one thing I would expect the boundary on the left to be 15 (x[0]=0), the right to be 5(x[0]=1). It looks as if the bottom right corner has a value of 5, i.e. as if the BC is rotated by -45 degrees. x[0]=1 and x[1]> 0.5 seems to have a value of 15 (should be 5)

It is 15 on the left and 5 on the right, at least in my plot. Are you looking at the second plot? There is no bc on Sw.

Mikael

>

> This would be easier if I can attach a screen shot, but I cannot find a

> way to do it.

>

> --

> You received this question notification because you are an answer

> contact for CBC.PDESys.

Yaakoub El Khamra (yelkhamra) said : | #5 |

The first plot is very weird for me. The second plot is what I am referring to. The second plot should be for Pw, and that's where I apply the BC. I think.

Should I just put the plot results online and paste a link or is there a way to attach a file in launchpad?

Yaakoub El Khamra (yelkhamra) said : | #6 |

sorry, hit problem solved, again.

Just send it to <email address hidden>

Mikael

Den Oct 22, 2012 kl. 8:45 PM skrev Yaakoub El Khamra:

> Question #211982 on CBC.PDESys changed:

> https:/

>

> Status: Answered => Solved

>

> Yaakoub El Khamra confirmed that the question is solved:

>

> The first plot is very weird for me. The second plot is what I am

> referring to. The second plot should be for Pw, and that's where I apply

> the BC. I think.

>

> Should I just put the plot results online and paste a link or is there a

> way to attach a file in launchpad?

>

> --

> You received this question notification because you are an answer

> contact for CBC.PDESys.

Yaakoub El Khamra (yelkhamra) said : | #8 |

Sent. Thank you very much in advance.

Yaakoub El Khamra (yelkhamra) said : | #9 |

I just sent you an email sir.

using: DirichletBC(

Regards

Yaakoub