1: //
   2: // Heccer : a compartmental solver that implements efficient Crank-Nicolson
   3: // integration for neuronal models.
   4: //
   5: 
   6: //////////////////////////////////////////////////////////////////////////////
   7: //'
   8: //' Heccer : testbed C implementation
   9: //'
  10: //' Copyright (C) 2006-2008 Hugo Cornelis
  11: //'
  12: //' functional ideas .. Hugo Cornelis, hugo.cornelis@gmail.com
  13: //'
  14: //' coding ............ Hugo Cornelis, hugo.cornelis@gmail.com
  15: //'
  16: //////////////////////////////////////////////////////////////////////////////
  17: 
  18: 
  19: #include <stdlib.h>
  20: 
  21: #include "../../heccer/addressing.h"
  22: #include "../../heccer/compartment.h"
  23: #include "../../heccer/heccer.h"
  24: #include "../../heccer/mathcomponent.h"
  25: #include "../../heccer/output.h"
  26: 
  27: 
  28: /* #define HECCER_TEST_INTERVAL_DEFAULT_START                   (-0.1) */
  29: /* #define HECCER_TEST_INTERVAL_DEFAULT_END                     (0.05) */
  30: /* #define HECCER_TEST_INTERVAL_DEFAULT_ENTRIES                 3000 */
  31: #define HECCER_TEST_REPORTING_GRANULARITY 10000
  32: #define HECCER_TEST_STEPS 40000
  33: #define HECCER_TEST_TESTED_THINGS ( HECCER_DUMP_VM_COMPARTMENT_MATRIX \
  34:                                     | HECCER_DUMP_VM_COMPARTMENT_DATA \
  35:                                     | HECCER_DUMP_VM_COMPARTMENT_OPERATIONS \
  36:                                     | HECCER_DUMP_VM_MECHANISM_DATA \
  37:                                     | HECCER_DUMP_VM_MECHANISM_OPERATIONS \
  38:                                     | HECCER_DUMP_VM_SUMMARY \
  39:         )
  40: #define HECCER_TEST_TIME_STEP (5e-7)
  41: 
  42: 
  43: struct Compartment compSoma =
  44: {
  45:     //m administrative overhead
  46: 
  47:     {
  48:         //m type of structure
  49: 
  50:         MATH_TYPE_Compartment,
  51:     },
  52: 
  53:     //m index of parent compartment, -1 for none
  54: 
  55:     -1,
  56: 
  57: /*     //m first mechanism */
  58: 
  59: /*     NULL, */
  60: 
  61: /*     //m number of mechanisms */
  62: 
  63: /*     0, */
  64: 
  65:     //m descriptive values, alphabetical order
  66: 
  67: /*     double dCm; */
  68: 
  69:     2.827430998e-11,
  70: 
  71: /*     double dEm; */
  72: 
  73:     -0.05939999968,
  74: 
  75: /*     double dInitVm; */
  76: 
  77:     -0.070,
  78: 
  79: /*     double dInject;           */
  80: 
  81:     3e-10,
  82: 
  83: /*     double dRa; */
  84: 
  85:     12732.40625,
  86: 
  87: /*     double dRm; */
  88: 
  89:     117891472,
  90: };
  91: 
  92: 
  93: //v a simple fast sodium channel
  94: 
  95: struct ChannelActInact caiNa =
  96: {
  97:     //m administrative overhead
  98: 
  99:     {
 100:         //m type of structure
 101: 
 102:         MATH_TYPE_ChannelActInact,
 103:     },
 104: 
 105:     //m first set of descriptive values, alphabetical order
 106: 
 107:     //m initial reversal potential
 108: 
 109:     0.045,
 110: 
 111:     //m get reversal potential from this intermediary, -1 for none
 112: 
 113:     -1,
 114: 
 115:     //m maximal conductance when all channels are permissive
 116: 
 117:     3.392917279e-06,
 118: 
 119:     //m contributes to this concentration pool, -1 for none, boolean indicator only.
 120: 
 121:     -1,
 122: 
 123:     //m activation description
 124: 
 125:     {
 126:         //m power, for a standard heccer, something between 1 and 4 or so.
 127: 
 128:         3,
 129: 
 130:         //m gate definition
 131: 
 132:         {
 133:             //m initial value, commonly forward over backward steady states
 134: 
 135:             0.05293248736,
 136: 
 137:             //m corresponding index in tables, set to -1 for initialization.
 138: 
 139:             -1,
 140: 
 141:             {
 142:                 //m forward kinetiks, commonly denoted with alpha or non-perm to perm rate
 143: 
 144:                 {
 145:                     //m multiplier
 146: 
 147:                     -4500,
 148: 
 149:                     //m multiplier membrane dependence, 0.0 for no dependence
 150: 
 151:                     0.1e6,
 152: 
 153:                     //m 2b: multiplier membrane dependence offset, 0.0 for no dependence
 154: 
 155:                     0.0,
 156: 
 157:                     //m choose between nominator or denominator, 1 means nominator, -1
 158:                     //m means denominator
 159: 
 160:                     -1.0,
 161: 
 162:                     //m nominator or denominator offset
 163: 
 164:                     -1.0,
 165: 
 166:                     //m membrane offset
 167: 
 168:                     0.045,
 169: 
 170:                     //m denormalized time constant
 171: 
 172:                     -0.01,
 173:                 },
 174: 
 175:                 //m backward kinetiks, commonly denoted with beta or perm to non-perm rate
 176: 
 177:                 {
 178:                     //m multiplier
 179: 
 180:                     4e3,
 181: 
 182:                     //m multiplier membrane dependence, 0.0 for no dependence
 183: 
 184:                     0.0,
 185: 
 186:                     //m 2b: multiplier membrane dependence offset, 0.0 for no dependence
 187: 
 188:                     1.0,
 189: 
 190:                     //m choose between nominator or denominator, 1 means nominator, -1
 191:                     //m means denominator
 192: 
 193:                     -1.0,
 194: 
 195:                     //m nominator or denominator offset
 196: 
 197:                     0.0,
 198: 
 199:                     //m membrane offset
 200: 
 201:                     0.07,
 202: 
 203:                     //m denormalized time constant
 204: 
 205:                     18e-3,
 206:                 },
 207:             },
 208:         },
 209:     },
 210: 
 211:     //m inactivation description
 212: 
 213:     {
 214:         //m power, for a standard heccer, something between 1 and 4 or so.
 215: 
 216:         1,
 217: 
 218:         //m gate definition
 219: 
 220:         {
 221:             //m initial value, commonly forward over backward steady states
 222: 
 223:             0.5961207962,
 224: 
 225:             //m corresponding index in tables, set to -1 for initialization.
 226: 
 227:             -1,
 228: 
 229:             {
 230:                 //m forward kinetiks, commonly denoted with alpha or non-perm to perm rate
 231: 
 232:                 {
 233:                     //m multiplier
 234: 
 235:                     70,
 236: 
 237:                     //m multiplier membrane dependence, 0.0 for no dependence
 238: 
 239:                     0.0,
 240: 
 241:                     //m 2b: multiplier membrane dependence offset, 0.0 for no dependence
 242: 
 243:                     0.0,
 244: 
 245:                     //m choose between nominator or denominator, 1 means nominator, -1
 246:                     //m means denominator
 247: 
 248:                     -1.0,
 249: 
 250:                     //m nominator or denominator offset
 251: 
 252:                     0.0,
 253: 
 254:                     //m membrane offset
 255: 
 256:                     0.07,
 257: 
 258:                     //m denormalized time constant
 259: 
 260:                     0.02,
 261:                 },
 262: 
 263:                 //m backward kinetiks, commonly denoted with beta or perm to non-perm rate
 264: 
 265:                 {
 266:                     //m multiplier
 267: 
 268:                     1.0e3,
 269: 
 270:                     //m multiplier membrane dependence, 0.0 for no dependence
 271: 
 272:                     0.0,
 273: 
 274:                     //m 2b: multiplier membrane dependence offset, 0.0 for no dependence
 275: 
 276:                     1.0,
 277: 
 278:                     //m choose between nominator or denominator, 1 means nominator, -1
 279:                     //m means denominator
 280: 
 281:                     -1.0,
 282: 
 283:                     //m nominator or denominator offset
 284: 
 285:                     1.0,
 286: 
 287:                     //m membrane offset
 288: 
 289:                     0.04,
 290: 
 291:                     //m denormalized time constant
 292: 
 293:                     -10.0e-3,
 294:                 },
 295:             },
 296:         },
 297:     },
 298: };
 299: 
 300: 
 301: struct ChannelAct caK =
 302: {
 303:     //m administrative overhead
 304: 
 305:     {
 306:         //m type of structure
 307: 
 308:         MATH_TYPE_ChannelAct,
 309:     },
 310: 
 311:     //m first set of descriptive values, alphabetical order
 312: 
 313:     //m initial reversal potential
 314: 
 315:     -0.08200000226,
 316: 
 317:     //m get reversal potential from this intermediary, -1 for none
 318: 
 319:     -1,
 320: 
 321:     //m maximal conductance when all channels are permissive
 322: 
 323:     1.017875206e-06,
 324: 
 325:     //m contributes to this concentration pool, -1 for none, boolean indicator only.
 326: 
 327:     -1,
 328: 
 329:     //m activation description
 330: 
 331:     {
 332:         //m power, for a standard heccer, something between 1 and 4 or so.
 333: 
 334:         4,
 335: 
 336:         //m gate definition
 337: 
 338:         {
 339:             //m initial value, commonly forward over backward steady states
 340: 
 341:             0.3176768857,
 342: 
 343:             //m corresponding index in tables, set to -1 for initialization.
 344: 
 345:             -1,
 346: 
 347:             {
 348:                 //m forward kinetiks, commonly denoted with alpha or non-perm to perm rate
 349: 
 350:                 {
 351:                     //m multiplier
 352: 
 353:                     -600,
 354: 
 355:                     //m multiplier membrane dependence, 0.0 for no dependence
 356: 
 357:                     10000,
 358: 
 359:                     //m 2b: multiplier membrane dependence offset, 0.0 for no dependence
 360: 
 361:                     0.0,
 362: 
 363:                     //m choose between nominator or denominator, 1 means nominator, -1
 364:                     //m means denominator
 365: 
 366:                     -1.0,
 367: 
 368:                     //m nominator or denominator offset
 369: 
 370:                     -1.0,
 371: 
 372:                     //m membrane offset
 373: 
 374:                     60e-3,
 375: 
 376:                     //m denormalized time constant
 377: 
 378:                     -10.0e-3,
 379:                 },
 380: 
 381:                 //m backward kinetiks, commonly denoted with beta or perm to non-perm rate
 382: 
 383:                 {
 384:                     //m multiplier
 385: 
 386:                     125,
 387: 
 388:                     //m multiplier membrane dependence, 0.0 for no dependence
 389: 
 390:                     0.0,
 391: 
 392:                     //m 2b: multiplier membrane dependence offset, 0.0 for no dependence
 393: 
 394:                     1.0,
 395: 
 396:                     //m choose between nominator or denominator, 1 means nominator, -1
 397:                     //m means denominator
 398: 
 399:                     -1.0,
 400: 
 401:                     //m nominator or denominator offset
 402: 
 403:                     0.0,
 404: 
 405:                     //m membrane offset
 406: 
 407:                     70e-3,
 408: 
 409:                     //m denormalized time constant
 410: 
 411:                     80e-3,
 412:                 },
 413:             },
 414:         },
 415:     },
 416: };
 417: 
 418: 
 419: int piC2m[] =
 420: {
 421:     2,
 422:     -1,
 423: };
 424: 
 425: 
 426: struct MathComponentArray mca =
 427: {
 428:     //m number of math components
 429: 
 430:     2,
 431: 
 432:     //m math component data
 433: 
 434:     NULL,
 435: 
 436:     //m math component index, initialize to NULL
 437: 
 438:     NULL,
 439: 
 440: };
 441: 
 442: 
 443: struct Intermediary inter =
 444: {
 445:     //m compartment array
 446: 
 447:     1,
 448: 
 449:     &compSoma,
 450: 
 451:     //m all other mathematical components
 452: 
 453:     &mca,
 454: 
 455:     //m compartment 2 first mechanism number
 456: 
 457:     piC2m,
 458: };
 459: 
 460: 
 461: struct OutputGenerator * pogVm = NULL;
 462: 
 463: double *pdVm = NULL;
 464: 
 465: char pcStep[100] = "";
 466: 
 467: #include "main.h"
 468: 
 469: 
 470: int main(int argc, char *argv[])
 471: {
 472:     //- set default result : ok
 473: 
 474:     int iResult = EXIT_SUCCESS;
 475: 
 476:     //- determine intermediary size, and allocate
 477: 
 478:     struct MathComponentInfo *pmciNa = MathComponentInfoLookup(caiNa.mc.iType);
 479: 
 480:     struct MathComponentInfo *pmciK = MathComponentInfoLookup(caK.mc.iType);
 481: 
 482:     int iChars = pmciNa->iChars + pmciK->iChars;;
 483: 
 484:     void *pmc = calloc(sizeof(char), iChars);
 485: 
 486:     //- prepare the mechanism intermediary
 487: 
 488:     struct ChannelActInact *pcaiNa = (struct ChannelActInact *)pmc;
 489: 
 490:     *pcaiNa = caiNa;
 491: 
 492:     struct ChannelAct *pcaK = (struct ChannelAct *)&((char *)pcaiNa)[pmciNa->iChars];
 493: 
 494:     *pcaK = caK;
 495: 
 496:     //- link the intermediary
 497: 
 498:     mca.pmc = pmc;
 499: 
 500:     //- create output elements
 501: 
 502:     pogVm = OutputGeneratorNew("/tmp/output");
 503: 
 504: //d generate output of membrane potential each step
 505: 
 506: #define HECCER_TEST_INITIATE \
 507:     pdVm = HeccerAddressCompartmentVariable(pheccer, 0, "Vm"); \
 508:     OutputGeneratorAddVariable(pogVm, "Vm", pdVm)
 509: 
 510: //d generate output of membrane potential each step
 511: 
 512: #define HECCER_TEST_OUTPUT OutputGeneratorAnnotatedStep(pogVm, sprintf(pcStep, "%i", i) ? pcStep : "sprintf() failed")
 513: 
 514:     //- do the simulation
 515: 
 516:     simulate(argc,argv);
 517: 
 518:     //- finish the simulation output
 519: 
 520:     OutputGeneratorFinish(pogVm);
 521: 
 522:     //- add the simulation output to the program output
 523: 
 524:     WriteOutput("/tmp/output");
 525: 
 526:     //- return result
 527: 
 528:     return(iResult);
 529: }
 530: 
 531: 
 532: #define main(argc,argv) simulate(argc,argv)
 533: 
 534: //t this prototype can give warning and perhaps errors.
 535: 
 536: int main(int argc, char *argv[]);
 537: 
 538: 
 539: #include "main.c"
 540: 
 541: 
 542: