~logan/ubuntu/quantal/cadabra/fix-for-1030273

« back to all changes in this revision

Viewing changes to src/modules/substitute.cc

  • Committer: Bazaar Package Importer
  • Author(s): Gürkan Sengün
  • Date: 2010-02-03 10:30:30 UTC
  • mfrom: (1.1.1 upstream) (0.1.3 sid)
  • Revision ID: james.westby@ubuntu.com-20100203103030-euf26dyybbsl1okf
* New upstream version.
* Bump standards version.

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
1
/* 
2
2
 
3
3
        Cadabra: a field-theory motivated computer algebra system.
4
 
        Copyright (C) 2001-2009  Kasper Peeters <kasper.peeters@aei.mpg.de>
 
4
        Copyright (C) 2001-2010  Kasper Peeters <kasper.peeters@aei.mpg.de>
5
5
 
6
6
   This program is free software: you can redistribute it and/or
7
7
   modify it under the terms of the GNU General Public License as
209
209
        std::vector<iterator> subtree_insertion_points;
210
210
        while(it!=repl.end()) { 
211
211
                bool is_stripped=false;
 
212
//              tr.print_recursive_treeform(std::cerr, repl.begin());
212
213
 
213
214
//              For some reason 'a?' is not found!?! Well, that's presumably because _{a?} does not
214
215
//      match ^{a?}. (though this does match when we write 'i' instead of a?. 
222
223
                         }
223
224
 
224
225
                if(loc!=comparator.replacement_map.end()) { // name wildcards
225
 
//                      std::cerr << "rule : " << *((*loc).first.begin()->name) << " -> " << std::endl;
226
 
//                      std::cerr << "going to replace " << *it->name << " with " << *((*loc).second.begin()->name) << std::endl;
 
226
//              std::cerr << "rule : " << *((*loc).first.begin()->name) << " -> " 
 
227
//                                       << *((*loc).second.begin()->name) << std::endl;
 
228
//              std::cerr << "going to replace " << *it->name << " with " << *((*loc).second.begin()->name) << std::endl;
227
229
 
228
230
                        // When a replacement is made here, and the index is actually
229
231
                        // a dummy in the replacement, we screw up the ind_dummy
236
238
                        ind_dummy.erase(exptree(it));
237
239
 
238
240
                        str_node::bracket_t remember_br=it->fl.bracket;
239
 
                        if(is_stripped)
240
 
                                 it->name=(*loc).second.begin()->name;
 
241
                        if(is_stripped || (it->is_name_wildcard() && !it->is_index()) ) { 
 
242
            // a?_{i j k} type patterns should only replace the head
 
243
                                // TODO: should we replace brackets here too?
 
244
                                it->name=(*loc).second.begin()->name;
 
245
                                it->multiplier=(*loc).second.begin()->multiplier;
 
246
                                it->fl=(*loc).second.begin()->fl;
 
247
                                }
241
248
                        else
242
249
                                 it=tr.replace_index(it, (*loc).second.begin());
243
250
                        it->fl.bracket=remember_br;
264
271
                        }
265
272
                else ++it;
266
273
                }
 
274
//      tr.print_recursive_treeform(std::cerr, repl.begin());
267
275
 
268
276
        // If the replacement contains dummies, avoid clashes introduced when
269
277
        // free indices in the replacement (induced from the original expression)
385
393
        }
386
394
 
387
395
 
 
396
/*  bug: cadabra-34.
 
397
 
 
398
    Vary should take into account the depth of an object in a more clever way than
 
399
         is currently done. Consider an expression
 
400
 
 
401
         A \partial{ A C + B } + D A;
 
402
 
 
403
    and A->a, B->b etc. This should vary to
 
404
 
 
405
         a  \partial{ A C + B } + A \partial{ a C + A c + b } + d A + a D;
 
406
 
 
407
    Right now it produces a total mess for the partial derivative, because it does
 
408
         not understand that the depth counting for factors inside the partial involves
 
409
         knowing about the top-level product.
 
410
 
 
411
    So what we should do is introduce a 'factor depth' and a 'sum index', which equal
 
412
 
 
413
         A \partial{ A C + B } + D A;
 
414
         1           1 1   1     1 1
 
415
                        1           1 1   1     2 2
 
416
 
 
417
    For all 
 
418
 */
 
419
 
388
420
vary::vary(exptree& tr, iterator it)
389
421
        : algorithm(tr, it)
390
422
        {
398
430
bool vary::can_apply(iterator it) 
399
431
        {
400
432
        if(*it->name=="\\prod") return true;
 
433
        if(*it->name=="\\sum") return true;
401
434
        if(is_single_term(it)) return true;
 
435
        if(is_nonprod_factor_in_prod(it)) return true;
 
436
        const Derivative *der = properties::get<Derivative>(it);
 
437
        if(der) return true;
 
438
        der = properties::get<Derivative>(tr.parent(it));
 
439
        if(der) return true;
 
440
        const Accent *acc = properties::get<Accent>(it);
 
441
        if(acc) return true;
 
442
        acc = properties::get<Accent>(tr.parent(it));
 
443
        if(acc) return true;
402
444
        return false;
403
445
        }
404
446
 
405
447
algorithm::result_t vary::apply(iterator& it)
406
448
        {
407
 
        if(is_single_term(it)) { // easy: just vary this term
 
449
        const Derivative *der = properties::get<Derivative>(it);
 
450
        const Accent     *acc = properties::get<Accent>(it);
 
451
        if(der || acc) {
 
452
                vary vry(tr, this_command);
 
453
 
 
454
                sibling_iterator sib=tr.begin(it);
 
455
                while(sib!=tr.end(it)) {
 
456
                        iterator app=sib;
 
457
                        ++sib;
 
458
                        if(app->is_index()) continue;
 
459
                        if(vry.can_apply(app)) {
 
460
                                vry.apply(app);
 
461
                                }// not complete: should remove object when zero
 
462
                        }
 
463
                expression_modified=true;
 
464
                return l_applied;               
 
465
                }
 
466
 
 
467
        if(*it->name=="\\prod") {
 
468
                exptree result;
 
469
                result.set_head(str_node("\\expression"));
 
470
                iterator newsum=result.append_child(result.begin(), str_node("\\sum"));
 
471
                
 
472
                // Iterate over all factors, attempting a substitute. If this
 
473
                // succeeds, copy the term to the "result" tree. Then restore the
 
474
                // original. We have to do the substitute on the original tree so
 
475
                // that index relabelling takes into account the rest of the tree.
 
476
                
 
477
                exptree prodcopy(it); // keep a copy to restore after each substitute
 
478
                
 
479
                vary subs(tr, this_command);
 
480
                int pos=0;
 
481
                for(;;) { 
 
482
                        sibling_iterator fcit=tr.begin(it);
 
483
                        fcit+=pos;
 
484
                        if(fcit==tr.end(it)) break;
 
485
                        
 
486
                        iterator fcit2(fcit);
 
487
                        if(subs.can_apply(fcit2)) {
 
488
                                subs.apply(fcit2);
 
489
                                expression_modified=true;
 
490
                                
 
491
                                if(fcit2->is_zero()==false)
 
492
                                        iterator newterm=result.append_child(newsum, it);
 
493
 
 
494
                                // restore original
 
495
                                it=tr.replace(it, prodcopy.begin());
 
496
                                }
 
497
                        ++pos;
 
498
                        }
 
499
                if(expression_modified && tr.number_of_children(newsum)>0) {
 
500
                        it=tr.move_ontop(it, newsum);
 
501
                        cleanup_nests(tr, it);
 
502
                        cleanup_expression(tr, it);
 
503
                        }
 
504
                else { // varying any of the factors produces nothing, variation is zero
 
505
                        zero(it->multiplier);
 
506
                        expression_modified=true;
 
507
                        }
 
508
                return l_applied;
 
509
                }
 
510
 
 
511
        if(*it->name=="\\sum") { // call vary on every term
 
512
                vary vry(tr, this_command);
 
513
 
 
514
                sibling_iterator sib=tr.begin(it);
 
515
                while(sib!=tr.end(it)) {
 
516
                        iterator app=sib;
 
517
                        ++sib;
 
518
                        if(vry.can_apply(app)) {
 
519
                                vry.apply(app);
 
520
                                if(app->is_zero()) {
 
521
                                        tr.erase(app);
 
522
                                        }
 
523
                                }
 
524
                        else {
 
525
                                // remove this term
 
526
                                tr.erase(app);
 
527
                                }
 
528
                        }
 
529
                expression_modified=true;
 
530
                return l_applied;
 
531
                }
 
532
        
 
533
        der = properties::get<Derivative>(tr.parent(it));
 
534
        acc = properties::get<Accent>(tr.parent(it));
 
535
 
 
536
        if(der || acc || is_single_term(it)) { // easy: just vary this term by substitution
408
537
                substitute subs(tr, this_command);
409
538
                if(subs.can_apply(it)) {
410
539
                        if(subs.apply(it)==l_applied) {
414
543
                        }
415
544
                return l_no_action;
416
545
                }
417
 
 
418
 
        exptree result;
419
 
        result.set_head(str_node("\\expression"));
420
 
        iterator newsum=result.append_child(result.begin(), str_node("\\sum"));
421
 
 
422
 
        // Iterate over all factors, attempting a substitute. If this
423
 
        // succeeds, copy the term to the "result" tree. Then restore the
424
 
        // original. We have to do the substitute on the original tree so
425
 
        // that index relabelling takes into account the rest of the tree.
426
 
 
427
 
        exptree prodcopy(it); // keep a copy to restore after each substitute
428
 
 
429
 
        substitute subs(tr, this_command);
430
 
        int pos=0;
431
 
        for(;;) { 
432
 
                sibling_iterator fcit=tr.begin(it);
433
 
                fcit+=pos;
434
 
                if(fcit==tr.end(it)) break;
435
 
 
436
 
                iterator fcit2(fcit);
437
 
                if(subs.apply_recursive(fcit2, false)) {
438
 
                        expression_modified=true;
439
 
                        iterator newterm=result.append_child(newsum, it);
440
 
                        // restore original
441
 
                        it=tr.replace(it, prodcopy.begin());
 
546
        
 
547
        if(is_nonprod_factor_in_prod(it)) {
 
548
                substitute subs(tr, this_command);
 
549
                if(subs.can_apply(it)) {
 
550
                        if(subs.apply(it)==l_applied) {
 
551
                                expression_modified=true;
 
552
                                return l_applied;
 
553
                                }
442
554
                        }
443
 
                ++pos;
444
 
                }
445
 
        if(expression_modified) {
446
 
                it=tr.move_ontop(it, newsum);
447
 
                cleanup_nests(tr, it);
448
 
                cleanup_expression(tr, it);
449
 
                }
450
 
        else { // substitute didn't act anywhere, variation is zero
 
555
                // else set to zero
451
556
                zero(it->multiplier);
452
557
                expression_modified=true;
 
558
                return l_applied;
453
559
                }
454
560
 
455
 
        return l_applied;
 
561
        return l_no_action;
456
562
        }
457
563
 
458
564