I am doing some numerical modeling for diffusion of a chemical through the
skin, and I use a for loop to time step my way through. For each loop I
solve a matrix and values in the Conc matrix are reset. Each present
solution depends on the solution from before, so I see no way to avoid a
for loop. I have this same function programmed on Matlab and it is so much
faster that the same function reprogrammed in S-PLUS. Does anyone have any
advice for making my for loop go faster. I am using S-PLUS 6 on Windows ME
and Windows XP. Below I have attached the main part of my program with the
for loop. Any help would be greatly appreciated. Thanks Alesia
for(t in dt:dt:ExposTime) {
BMatrix[1, 1] <- Conc[1, 1] *
((1/dt) - (A * K1 *
P1)/(2 * V)) + Conc[
2, 1] * ((A * K1)/
(2 * V))
BMatrix[2, 1] <- Conc[1, 1] *
((K1 * P1)/(2 * (dx1))) +
Conc[2, 1] * ((1/dt) -
(K1/2 * dx1) - (D1/
(2 * (dx1^2)))) + Conc[
3, 1] * (D1/(2 * (
dx1^2)))
BMatrix[3:(CSize1 - 1), 1] <-
(1 - 2 * Alpha1) * Conc[
3:(CSize1 - 1), 1] +
Alpha1 * Conc[4:CSize1,
1] + Alpha1 * Conc[
2:(CSize1 - 2), 1]
BMatrix[CSize1, 1] <- Conc[
CSize1, 1] * ((1/dt) -
(D1/(2 * (dx1^2))) -
(K2/(2 * (dx1)))) +
Conc[(CSize1 - 1),
1] * (D1/(2 * (dx1^
2))) + Conc[(CSize1 +
1), 1] * ((K2 * P2)/
(2 * dx1))
BMatrix[(CSize1 + 1), 1] <-
Conc[(CSize1 + 1),
1] * ((1/dt) - ((K2 *
P2)/(2 * dx2)) - (
D2/(2 * (dx2^2)))) +
Conc[(CSize1 + 2),
1] * (D2/(2 * (dx2^
2))) + Conc[CSize1,
1] * (K2/(2 * dx2))
BMatrix[(CSize1 + 2):(CSize -
1), 1] <- (1 - 2 *
Alpha2) * Conc[(CSize1 +
2):(CSize - 1), 1] +
Alpha2 * Conc[(CSize1 +
3):(CSize), 1] + Alpha2 *
Conc[(CSize1 + 1):
(CSize - 2), 1]
BMatrix[CSize, 1] <- Conc[
CSize, 1] * ((1/dt) -
(D2/(dx2^2))) + Conc[
(CSize - 1), 1] * (
D2/(2 * (dx2^2)))
#Now we need to get the solution for the ConcVector, AMatrix
Defined Above
Conc <- solve(AMatrix, BMatrix)
#
SkinMassPerArea <- ((sum(Conc[
2:CSize1]) * dx1) +
(sum(Conc[(CSize1 +
1):CSize]) * dx2))
MassAbsorbed <-
MassLostToBloodStream +
(A * SkinMassPerArea)
SCMASS <- A * (sum(Conc[2:
CSize1]) * dx1)
VEMASS <- A * (sum(Conc[(CSize1 +
1):CSize]) * dx2)
PercentInSC <- (SCMASS/
OriginalVehicleMass) *
100
PercentInVE <- (VEMASS/
OriginalVehicleMass) *
100
PercentInBlood <- (
MassLostToBloodStream/
OriginalVehicleMass) *
100
MassLostToBloodStream <-
MassLostToBloodStream +
(A * dt * D2 * (Conc[
CSize, 1]/dx2))
MassBalance <- (Conc[1, 1] *
V) + MassAbsorbed
}
|