dsb-lab / CellBasedModels.jl

Julia package for multicellular modeling
MIT License
18 stars 3 forks source link

Questions on medium secretion/Neumann boundary conditions #34

Closed dirypan closed 11 months ago

dirypan commented 11 months ago

In the tutorial of the documentation, Defining an ABM model has a diffusive medium example

model = ABM(2,
    model = Dict(
        :D => Float64
        ),
    agent = Dict(
        :secrete => Float64
    ),
    medium = Dict(
        :p=> Float64  #Number of neighbors that the agent has
        ),
    agentRule = quote
        p += secrete/dt*dx*dy #Add to the medium the secretion content of each neighbor at each time point
    end,
    mediumODE = quote
        if @mediumInside()
            dt(p) = @∂2(1,p)        #Diffusion inside
        elseif @mediumBorder(1,-1)
            p = p[2]                #Newman (reflective) boundaries on the borders
        elseif @mediumBorder(1,1)
            p = 0                   #Dirichlet (absorvant) boundary
        end
    end

I have two questions related to the example since I also want to implement a diffusive system.

  1. In the agent rule, I think the correct rule should be
    p += secrete*dt/(dx*dy)

    where "secrete" is the secretion rate of one agent. This is consistent with units and rescaling of dt,dx,dy.

  2. I don't understand the reflective boundary conditions
    p = p[2]                #Newman (reflective) boundaries on the borders

    This doesn't make sense to me, and if I run this code with evolve it will raise an error. I do want to implement a reflective boundary condition, but I don't know how to access neighboring medium grid points so I can manually write down the condition like here

gatocor commented 11 months ago

Hi, thanks for trying the package!

  1. You are right. Thanks for pointing it out.
  2. The Newman boundary conditions I was considering are of the type $$\partialx f(x) = 0$$ which discretising (for example in forward pass), $$(f{i+1}-f{i})/\Delta t = 0$$ and reorganizing we have, $$f{i+1}= f_{i}$$ which will have to be evaluated at boundary conditions only.

For accessing any point in the medium, you can always access as an array the medium variable. The correct way of accessing it should be:

    mediumODE = quote
        if @mediumInside()
            dt(p) = @∂2(1,p)        #Diffusion inside
        elseif @mediumBorder(1,-1)
            p = p[2,j]                #Newman (reflective) boundaries on the borders
        elseif @mediumBorder(1,1)
            p = 0                   #Dirichlet (absorvant) boundary
        end
    end

However, I am seeing there is a bug on the code as it overrides the '[2,i]' indexer. I am correcting it now. I will let you know when it is corrected (hopefully today or tomorrow should be done).

gatocor commented 11 months ago

I have solved the bug and updated the documentation @dirypan

I have made a new version of the code that solves that and other problem. Both are tagged in the new release of the program.

Let me know if this have solved your problems.

dirypan commented 11 months ago

Thank you @gatocor ! I think it solves the problem, now I can access the neighboring medium array by

p = p[2,i2_]

Thank you again for developing this great tool. Here is the experiment and model I am trying to simulate with the tool, I will refer to this package if I successfully build some new understanding of the model.